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FOREWORD 
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LL90A1760K with Lockheed Missiles & Space Company, Inc.; funded by the 
Langley Research Center of the National Aeronautics and Space Administration. 
The Contracting Officer’s Technical Representative was J. C, Robinson. 

The purpose of this document is to describe the internal logic of 
designated areas of SPAR system code. In anticipation of future additions, 
the following organization has been chosen; the document will consist of a 
sequence of named chapters, with each chapter incorporating its own table of 
contents and page/section notation. Chapters are tabulated in a General 
Table of Contents. 
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GENERAL INFORMATION 


Readers of this section should have a thorough understanding of the assumed 

6 - 

stress field-minimum complementary energy method for computing finite element 

stiffness and stress matrices. This method was originally proposed in 1964 

* 

by T-H.H. Plan, Reference 1 . The basic theoretical procedure and the assump- 
tions associated with the SPAR/K2D implementation of the Plan method are 

* 

presented in Reference 2 . This section documents the K2D-series of subroutines 
for computing stiffness and stress matrices for two-dimensional finite elements. 

K2D can accommodate the following: 

• 3- and 4-node element geometries, 

• membrane, plate bending, uncoupled membrane + bending, and 
coupled membrane + bending constitutive relations, and 

e quadrilateral elements with node 4 positioned slightly outside 
the x-y plane of the element reference frame. 

As discussed in References 1 and 2, the following three principal computa- - 
tional tasks are required: 

1. Compute the Internal strain energy matrix, H, such that 

V = H B, 

2. Compute the boundary work matrix, T, such that 

W = B* T q, and 

3. Form the element stress matrix = H and the element stiffness 

* 

matrix = T H T. 

* Ref. 1. Plan, "Derivation of Element Stiffness Matrices by 

Assumed Stress Distributions, AIAA J. 2, 1333-1336 (1964). 

Ref. 2. SPAR Reference Manual, Vol. IT, Section A. 
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The K2D-serxes of subroutines are organized as shown below 


Element Geometry 


Stress-Strain Relation 


I 

J 


K2D 




Assumed Stress Distribution 


Degree-of -Freedom Map 


1 


HGEN 

driver routine 

HGEND 

compute membrane 
and/or bending 
energies 

HMBGEN 

compute membrane- 
bending coupling 

ITQUAD ] 

n m 

xntegrate x y 

DSUM 

over element area 

SYMVRT 

invert H 


TGEN 

. driver routine 

TTGEN 


TT6X3 

set up boundary 

TT10X3 

> stress resultant 

RGEN 

TCB 

functions 

AGEN 

set up boundary 

ATD 

deform, functions 

NCALNA 

integrate and 
transform to 

CXA 

element coordinates 


Boundary Work Matrix, 

T 

1 

Node 4 lies out 
of the x-y plane 

i 


Flat 
Element 


WARPT 


Transformed Boundary Work, 
T — 


Strain Energy Matrix Inverse, 



SSTM 


' -1 

Stress Matrix = H T 
Stiffness Matrix = T*H T 
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K2D (MAIN ROUTINE) 


The following inputs are supplied through the K2D callxng sequence: 

A = Upper triangular portion (stored by column) of the 
element stress-strain relation. The following 
equations i-l-lustrate the required A configurations: 

Membrane Elements 
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Shear Panels 



a, e 
1 xy 


( 5 ) 


XX(2, NSIDES), XX(I,J) = directxon-I position of node J, 
relative to the element reference frame 

NSIDES, either 3 or 4 

NBM = number of assumed membrane stress resultant coefficients 
(maximum = 10) 

NBB = number of assumed bending stress resultant coefficients 
(maximxjm - 15) 

ISHELL = 0, no membrane-bending coupling 

0, membrane-bending coupling in A 

X34 = out-of-plane Z location of node 

1X34 = 0, no out-of -plane transformation 

> 0, transform T to account for non-zero X34 


The following arrays are returned through the K2D calling sequence: 

S = upper triangular portion (stored by column) of the NQT x NQT 
element stiffness matrix, 

STM (NBT.NQT) = element stress matrix, 

NBT = NBM + NBB, and 

NQT = order of the stiffness matrix. 


The following error conditions are returned through the common block ERRK2D 
(K2D line number 10) . 


AERR 

NERR 

Meaning 

4HISNG 

n 

H matrix singularity, row n 

4HNEGD 

n 

n negative diagonal terms 
encountered computing 

4HC0RE 

n 

n additional core locations 
required in array W 

Blank 

o 

no error 
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The procedure incorporated in K2D is outlined below. 


Step A. (K2D line numbers 20 through 31) 


Set up the assumed membrane stress resultant distribution 
via data statements. The arrays IBM, BM, IXYM, and XYM 
are used to characterize the equation below. 



For the preceding equation, IBM and BM characterize the 
matrix of b's as 


IBM = p3 6 4 10 8 0 *1 , and 

2 5 7 9 10 0 

[l 7 6 0 0 10 _ 

BM = ri. 1. 1. 1. 1. 0."J . 

1 . 1 . 1 . 1 . 1 . 0 . 

Ll. -i. -1. 0. 0. -2. 


XYM and IXYM characterize the vector of x’s and y’s. For 
the term, IXYM (1,J) is the exponent of x, IXYM (2,J) is 
the exponent of y, and XYM (J) is the coefficient. For 
example: 


Vector IXYM(1,J) IXYM(2,J) XYM(J) 
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If P(I) 

XYM(I) 

IXYM(1,J) 

IXYM(2,J) 


I — membrane stress, resultant 

(x.e. P(l) = N , P(2) = N P(3) = N ) , 

X y xy 

k, 

Z, 

n, and 
m, then 


P(I) = P(I) + b^ * C Z y“, for J == 1, 6. (7) 

For current applications, Z=1.0 for all values of J, This 
is evident from the dataed array XYM. Provisions for two 
constants, C and Z, are included for possible future 
applications. 

Current dimensions permit a maximum of 10 membrane stress 
resultant coefficients. The input variable NBM controls 
the selection of 10 or less. Note that 

NBM = 1 = pure shear stress distribution, used 
for E44 elements, and 

NBM = 3 = constant stress distribution, used for 
E31 elements. 


Step B . (K2D line numbers 33 through 48) 


Set up the assumed plate bending stress resultant distri- 
buting via data statements. Similar to the membrane 
representation, BB, IBB, XXB, and IXYB characterize the 
equation below. , ^ 

1 
X 


X 


M 

X 



0 

^6 

0 

^4 

0 

h4 

h2 

ho 

0 

M 

y 


^2 

0 


0 


0 

h3 

hs 

hi 

0 

M 

xy 

KS 


0 


0 

"8 

0 

0 

0 

h4 

h5 

Qx 



-^8 

^4 


ho 

0 

0 

0 

0 

0 

1 

o 

1 




^1 

0 

h5 

-h4 

0 

0 

0 

0 


y 

y 


y 

xy 

xy 


( 8 ) 
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Step C. 


(K2D line number 50) 


Set up degree~of-freedom (DOF) indicators 

NQM (1) = 3 = 3-node membrane DOF's 

NQM (2) = 5 = 4-node membrane DOF's 

NQB (1) = 6 = 3-node plate element DOF's 

NQB (2) = 9 ~ 4-node plate element DOF's 

t 

Step D. (K2D line numbers 52 through 61) 

Set up MAPM(4,5) = membrane DOF map and MAPB(4,5) 
== plate bending DOF map. 

MAPM(I,J) = membrane degree-of-freedom associated 
with direction-J motion of node I. 

MAPB(I,J) = bending degree-of-freedom associated 
with direction-J motion of node I. 

J = 1, ?, 3 = displacement in direction-J 

J = 4, 5 = rotation about axis J-3. 


Step E, (K2D line numbers 62 through 64) 

Initialize error indicators, AERR and NEKR. 

IWARP = 0 eliminates the possibility of activating 
alternative boundary work calculations concernxng 
warped boundaries. All warped element considerations 
are handled by the routine WARPT. 


Step F. (K2D line numbers 66 through 71) 

Set up the nodal intrinsic position coordinate array, 
X(2,4). For 3-node elements: 

X(I,4) = X(I,3), 1= 1,2. 

Step G. (K2D line numbers 72 through 99) 

Set up the following quantities: 

NBOTH =0, membrane or bending alone 
NBOTH = 1, membrane plus bending 
NBT = total number of assumed stress 
resultant coefficients 
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NXM 

NXB 

NQT 

MAPQ(4,5) 


order o£ the vector XYM 

order of the vector XfB 

number of element degrees of freedom 

element degree-of -freedom map 

MAPQ(I,J) = element DOF associated 

with directlon-J motion of node I. 


Step H. (K2D line numbers 107 and 108) 

Compute the element strain energy matrix inverse, H 
by calling HGEN. is stored in the array W 

beginning at location IH, 


Step I. (K2D line numbers 111 through 117) 


Set up NI(4) and NJ(4). Side K of the element is 
bounded by nodes NI(K) and NJ(K). The following table 
summarizes NI and NJ. 

3~Node Elements 4-Node Elements 

K NI(K) NJ(K) NI(K) NJ(K) 


113 

2 3 2 

3 2 1 

4 _ _ 


1 

4 

3 

2 


4 

3 

2 

1 


Step J. (K2D line number 120) 

Compute the element boundary work matrix, T, by calling 
TGEN. T is stored in the array W beginning at location 
IT. 


Step K. (K2D line niamber 123) 

Modify T to account for an out-of-plane node 4 by calling 
WAEPT. 


Step L. (K2D line number 126) 

Call SSTM to compute the intrinsic stress and stiffness 
matrices , 
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ELEMENT STRAIN ENERGY 


For a coupled membrane - plate bending fxnite element formulation, 
the stress-strain relation for the cross-section is written as 



The assumed stress resultant distributions for the element are 


F 

m 


P B , and 
m m 


P, B, , where 
b b’ 


B “ 
ra 


(bj^ b^ •• 


B, = 


(b^ ^2 •• 


vector of membrane stress 
resultant coefficients 

vector of bending stress 
resultant coefficients. 


( 10 ) 
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NBM and NBB represent the number of membrane and bending stress resultant 

coefficients assumed for the formulation. All elements of P and P, are 

m D 

of the general form x y™. 

The element strain energy is expressed in terms of the following area 
integral 



H , H , and H,, are submatrices of the desired strain energy matrix H, 
mm’ mb i>b 
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For a membrane-plate bending coup'led element, H is computed according 
to th€ following chart. 


From K2D 

Assumed Stress' Distributions, 
Stress-Strain Relations, 
Element Geometry 



SYMVRT 


H 


'-1 


For a pure membrane, plate bending, or uncoupled element formulation only 
applicable portions of the above chart are followed. 
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HGEN (PARENT ROUTINE = K2D) 


The following inputs are supplied through the HGEN calling sequence; 

NBM, NBB = number of assumed membrane and bending stress 
resultant coefficients 

X(2,4) « nodal position coordinates 

AA = upper triangular portion (stored by column) of 

the stress-strain matrix 

XYM, IXYM, XYB, IXYB, IBM, BM, IBB, BB = 

arrays characterizing the assiimed membrane and 
bending stress resultant distributions (see K2D) 

NXM, NXB = order of the vectors XYM and XYB 

ISHELL = 0, no membrane-bending coupling in AA 
0, membrane-bending coupling in AA 

NBT = NBM + NBB 

HGEN directs the computation of the element strain energy matrix inverse 
H(NBT,NBT). H is returned through the HGEN calling sequence. 

The procedure Incorporated in HGEN is outlined below. 

Step A. (HGEN line number 8) 

Set up via the data statement: 

NC =3, order of the Aj, and A^ 
stress-strain matrices 
NSM =3, number of rows in IBM and BM 
NSB =5, number of rows in IBB and BB 

Step B. (HGEN line numbers 9 and 10) 

Initialize the IQT array to 0. In the common 
block /CK204/ 
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QT(I,J) = 


(I-l) (J-1) 

X y dArea . 


Area 


If IQT(I,J) = 0, QT(I,J) is not computed. 

7 ^ 0, ,QT(I,J) is computed. 


Step Ct (HGEN line nxiraber 11 through 13) 

XX(I) and Y5f(I) are the x and y positrons 
relative to the element reference frame of 
node I, 


Step D. (HGEN l ine numbers 15 through 17) 
Initialize H to zero. 


Step E. (HGEN line numbers 19 through 49) 

Set up A(3,3,3), the three stress-strain matrices 
stored in rectangular fashion. Entries to A are 
extracted from the input array AA. See Eq. 9. 


El. (HGEN line nximber 25) 

For NBM=1 a shear panel formulation is implemented, 
and the stress-strain relation consists of a single 
term as shotm in the following equation: 



E2. (HGEN line 
A(I,J,1) = 


= A(l) e 

xy 

numbers 27 through 33) 
membrane stress-strain matrix 


E3. (HGEN line numbers 34 through 42) 

A(I,J,3) - bending stress-strain matrix 


E4. (HGEN line numbers 43 through 48) 

A(I,J,2) = coupling stress-strain matrix 
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Step F, (HGEN line numbers 53 through 59) 

Compute the membrane contribution to H by- 
calling HGEND. Locations within the temporary 
array W are assigned for storing 

Step G. (HGEN line numbers 61 through 67) 

Compute the bending contribution to H by 
calling HGEND. Locations within the temporary 
array W are assigned for storing P^. 

Step H. (HGEN line numbers 69 and 70) 

Compute the membrane-bending coupling contri- 
bution to H by calling HMBGEN. 

Step I« (HGEN line number 73) 

Invert the H matrix by calling SYMVRT. 
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HGEND (PAEIENT ROUTINE = HGKN) 


The following inputs are supplied through the HGEND callxng sequence; 


NB 


INCB 


NBT 


= number of assumed stress coefficients. 
For membrane energy, NB = NBM. 

For plate-bending energy, NB = NBB. 

= stress coefficient increment. 

For membrane energy, INCB = 0. 

For plate-bending energy, INCB = NBM. 

= total number of assumed stress 
Coefficients = NBM + NBB. 


C(NC,NC) = stress-strain'matrix 


XY(NX), IXY(NX) 1 define the assumed stress 
IB(NS,NX), B(NS,NX)J distribution (see K2D) 

INCC = increment for locating the beginning 
stress resultant row number within B 
and IB. INCC = 0 for current applications. 

The followl,ng arrays are returned through the HGEND calling sequence; 


IPX 

IPY 

P 


These arrays define the assumed stress 
distribution in a configuration that is 
convenient for computing H terms. 


H(NBT,NBT). Contributions to the element strain 
energy are added to this array. 


K2D-16 




The procedure incorporated in HGEND is outlined below: 


Step A. (HGEND line numbers 11 through 26) 


The input stress-resultant distribution 
matrices (XY, IXY, B, and IB) are 
rearranged for computing strain energy 
terms . The input matrices characterize 
the stress distribution as follows (see 
K2D discussion): 


s. 


1 




2 

• 


• 

s 


L NSj 



^i • * 


NB 


IB and B 


XY and IXY- 


(14) 


The following arrangement is used for 
computing R, Eq. 10. 


s. 


1 


s.. 

— 

2 

» 


• 

s 


L NSj 

_ 


P P 

11 12 

P P 

21 22 


L NS,1 


P, IPX, and IPY 



(15) 


In the above equation b^ is the a- 


th 


coefficient and P^ 

n ra 
a X y . 


stress 
IS of the general form 
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The terms of P, IPX, and IPY are extracted 
from B, IB, XY, and IXY. If = a x“ y“, 
then 

P(I,J) - a 
IPX(I,J) = n 
IPY(I,J) = m. 

Step B. (H6END line numbers 28 through 47) 

Compute contributions to H of the form 

P CP dArea, see Eqs. 11 and 13. 

Area 


Bl. (HGEND line numbers 29 through 31) 


IH and JH are the row/col\amn designation 
within H. 



B2. (HGEND line numbers 34 through 44) 

A typical term in the integral is 

H(IH,JH) = H(1H,JH) + PC y x^ y^ dArea 

Area 


(16) 
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If 



MY 

y 


dArea has not been previously 


computed, IIQUAD is called. The common 
block /CK2D1/ IS used to transfer the 
necessary data to and from ITQUAD. 
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HMBGEN (PARENT ROUTINE HGEN) 


The following inputs are supplied through the HMBGEN calling sequence: 


NBM and NBB = number of assumed membrane and bending 
stress-resultant coefficients 


C(NCM,NCB) = matrix relating membrane strains to 
bending stress resultants 


IPXM 

IPYM 

PM 


arrays characterizing the assumed 
membrane stress-resultant distribution, 
see HGEND 


IPXB 

IPYB 

PB 


arrays characterizing the assumed 
bending stress-resultant distribution 


NBT == NBM + "NBB = order of the H matrix 


Contributions to the element strain energy associated with coupled membrane- 
bending action are computed by HMBGEN and inserted in H(NBT,NBT). 


The procedure incorporated by HMBGEN is outlined below: 


Step A. (HMBGEN line numbers 12 through 30) 


Compute contributions to H of the 
form 


f* 

J ” 

Area 


P CP, dArea, see Eq. 12. 
m ^ ^ 
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Al. (HMBGEN line numbers 13 through 15) 

IH and JH are the row/column designation 
within H. 





IH-I 

JH=J+NBM 


A2. (HMBGEN line numbers 18 through 28) 
A typical term in the integral is 


H(IH,JH) = H(IH,JB) + PC 



MY 

y 


dArea 


As in HGEND, the value of the integral 
is taken from QX of the common block 
/CK2D4/ 


(17) 
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SYMVRT (PABEHT ROUTINE = HGEN) 


The following inputs are supplied through the SYMVRT callxng sequence: 

A(N,N) = S3nmnetric matrix 

NHD = dimension of A. The maximum 
for NHD is 30 due to current 
dimensions of the temporary 
arrays B and T. 

SYMVRT computes A ^ and stores the result in the space occupied by the 
original A matrix. 

Error conditions encountered during inversion are returned through the 

t 

common block /ERRK2D/. If AEEtR = ISNG, NERR contains the row nimiber of 
A at which a singularity is encountered. This constitutes an exit 
condition; the matrix inversion is terminated. If AERR = NEGD, NERR 
contains the number of negative diagonal terms encountered during 
inversion. This indicates that A is not positive definite. 
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QUADRILATERAL AREA INTEGRATOR 


A procedure is developed below for integrating a function of the general 
form y™ over an arbitrarily shaped quadrilateral surface. 


The quadrilateral area is illustrated below, 
y 



The quadrilateral Integral is computed as the summation three partial 
surface integrals . 




y™ dArea = A^ + + A^ 


(18) 


X, 0 
4 


m 


^3 -^4 


y dx dy ; y = + X - x ~ ^4^ 

3 4 


(19) 


^4 y 

A, - f j x" 


^ J J 

y ax dy ; y = — x 

X / 


0 0 


- 


( 20 ) 


X, 


2 y 


■//■■ 


m 


A 3 = / / x“ y “ dx dy 

^3 0 




( 21 ) 
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Aj^ and integrals are transformed to relative coordinates, so that 


where x = x - x. 


(x^-x^) 




// 


(x + x^)" y“ dx dy ; y = y^ 


y. - y 


x„ - X, 


4 - 

X 


( 22 ) 


where x = x - X 2 


A 


3 


(X3-X2) y 



(X 2 + x)^ y® dx dy 



(23) 


The general form of Aj^ and A^ is 


A == 


/ 


/ 


(Xj^ + x)“ y“ dx dy 


y = a + bx. 


(24) 


Integrating 


A = 


T 

J (Xj^ + x)“ (a + bx)”^^ dx. 
0 


(25) 
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In the preceding expression 


_ n _n 

(Xj^ + x) = dj X + X + ... and 


iirfl m+l tn 

(a + bx) = X + C2 X + ... c^2’ 


( 26 ) 


A is now rewritten as 


(n+nrf*l) 


n-hn 


A = 


nH-i 


djCjX 


+(djCj + 


d2C^)x 


+ . • 


■Vl 


(27) 


A2 and A^ are computed according to the above equation. is evaluated 
directly as 


n+i nri-1 

= (28) 
’2 (n+nH-2) (irH-l) 


The routines used to implement this quadrilateral integration procedure 
are ITQUAD and DSUM. 
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ITQUAB (PARENT ROUTINES = HGEND- and HMBGEN) 


The following- inputs are supplied through the labeled* common block 
/CK2D1/': 

X(4) I X(I) and Y(I) are the x and, y locations 
Y(4) j of node I 

N = exponent of x 

M - exponent of y 

ITQUAD directs the coinputation of 

AXNYM = dArea. 

Area 

AXNYM is returned in the calling routine' through /CK2D1/. 

The procedure incorporated in ITQUAD is outlined below:' 

Step A. (ITQUAD line numbers 7 and 8) 

Initialize AR(1) , AR(2) , and AR'(3) to zero*. 

Step B» (ITQUAD line nianbers 9 through 49) 

Integrate over zone 1*; compute AR(1) 

Bl. (ITQUAD line numbers 10 through 13) 

XL is a nominal element dimension against 
which other element dimensions are compared 
to determine whether or not they should be 
‘ considered zero. 

If D/XL< ZERO, D = 0. 

B2. (ITQUAD line number 14) 

If - x^ = 0, AR(1) - 0. This convention 
is used for integrating over triangular surfaces. 
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B3. (ITQUAD line numbers 15 and 16) 


Set up y = A + B X 


B4. (ITQUAD line numbers 20 through 22) 


For B = 0, 


B5. (ITQUAD line numbers 25 through 33) 

lH-1 _i«-l M 

For B > 0, y = C(l) x + C(2) x + ... C(M+2) 


B6. (ITQUAD line numbers 35 through 44) 

N N N-1 

(x^ + x) = CC(1) X + CC(2) X + . . CC(N+1) (29) 


B7. (ITQUAD line numbers 45 through 49) 


Compute the AR(1) integral by calling DSUM. 
Data IS transferred to and from DSUM through 
the common block /CK2D2/. 


SUM = 


/ r 1 r ^ 

C(l) X + ... CC(1) X + ... I 


dx 


(30) 
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Step C. (ITQUAD line numbers 51 through 54) 

Integrate over zone 2; compute AR(2) . 
If = 0, AR(2) = 0. 

Step D, (ITQUAD line numbers 56 through 75) 

Integrate over zone 3; compute AR(3) . 


Dl. (ITQUAD line numbers 59 through 68) 

N N N+1 

(x 2 + x) = CC(1) X + CC(2) X + ... CC(N+1) 


D2. (ITQUAD line numbers 69 through 73) 
Call DSUM to compute 


D3. 


T 

J- 


M+1 


SUM = J C(l) X 
0 



_N 

X + 


(ITQUAD line number 75) 
Compute AR(3) 


CC(N+2) dx 


Step E. (ITQUAD line numbers 77 through 81) 
Compute 

3 

1=1 


(31) 


(32) 


(33) 
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DSUM (PARENT ROUTINE = ITQUAD) 


The followxng inputs are supplied through the common block /CK2D2/ 


N 


} h = 


N 


N-1 


B(l) X + B(2) X + ... B(N+1) 


( 34 ) 


M 

MM 

C 


M 


M-1 


M+l-MM 


= C(l) X + C(2) X + ... C(MM) X 


(35) 


T = integrating limit 


DSUM evaluates the following integral and returns the result through 
/CK2D2/ . 

T 

SUM 


■ fh b 


0 


(36) 


substituting for f^ and and expanding the product yields 


SUM = 


N+M N+M-1 

B(l) C(l) X + B(l) C(2) + B(2) C(l) X 

r 1 _N+M-2 

+ B(l) C(3) + B(2) C(2) + B(3) C(l) x 
*- _M+1-MM] -* 

+ ... + B(Nhl) C(MM) X j dx 


(37) 


SUM is evaluated by taking each term of the above summation and 
integrating from x=0tox=T. A typical term of the sum is 


D = 


r 

u 

^ I=I 


IE 

B(I) * C(I) X 

=IB 


dx, where 


(38) 


J = N + M- l-IE+2, 1<J<MM. 


K2D-29 



BOUNDARY WORK 


The work performed by element stress resultants acting through boundary 
deformations is conq)uted as 


nsxdes 



where ¥n is the work performed along the n element boundary. 

Associated' with each element boundary is a boundary reference frame, 
X — y — Zj as shown below for boundary n, connecting nodes 1 and j . 



y is parallel with z; x extends from i through j, and z is directed 
toward the element interior. Points lying along i-j are located by 
X = t j + X COS0 and y = t^ + x sin0 . 

Functions characterizing the deformation of boundary n are 

u (x) = displacement in direction-z, 

V (x) = displacement in direction~x 
w (x) = displacement in direction-y, 

(j> (x) = rotation about axis z, and 
e (x) = rotation about axis x. 


K2D-30 


t 



The element stress resultants acting along boundary n of a coupled 
membrane plate-bending element are illustrated below. 



The boundary work performed along boundary n is 

L 

Wn = j" (N^ u + V + 0^ w + 0 0) dx. 

0 

where P, = N , P„ = N , P„ = M , P, = M , P^ = 0 > 

1 u 2 V 3 u’ 4 V 5 ^ 

= u, U 2 = V, TJ^ = (|), - 0j and = w, the above expression 

IS written as 

L 

5 r 

" Z / \ \ O*- <39) 

k=l J 
0 

For a pure membrane formulation k=l, 2. For pure plate-bending 
k=3, 4, 5. Por coupled membrane plus bending k=l through 5. 

Each boundary stress resultant, P^^, and each boundary deformation 
function, Uj^, is represented by a polynomial function in x. 
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as shown below 


\ - (G^ G2 ... Cj^^) 


X 

- 2 
X 


- NC-1 

X 


= C X , 
c ’ 


_NA-1 

Uj^ = { 1 X . . . X ) 


* 

X a. 
a 


(-40) 


(41) 


^NA 


The boundary work performed along boundary n may now be expressed 
as 



Where B represents the vector of undetermined stress -resultant 
coefficients and Q represents the vector of intrinsic nodal motions, 
the following transformations are performed; 


C = T , B and a = T Q. 
cb aq ^ 

(43) 

xp * * 

Wn=>BT,X T Q. 

/■ cb ca aq ^ 

k=i 
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Summing over all sides 


A 

W f= B TQ, where T = 


nsides 5 
n=l k=l 


T is the desired boundary work matrix. 



T is computed according to the following chart 


From K2D 

Assumed Stress Distribution 
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TGEN (PARENT ROUTINE = K2D) 


The followxng inputs are supplied through the TGEN calling sequence* 

IBM 1 element membrane stress-resultant distribution 

BM j matrices (see K2D) 

NBM = number of assumed membrane stress-resultant 
coefficients 

IBB 1 element plate-bending stress-resultant 

BB I distribution matrices (see K2D) 

NBB = number of assumed bending stpess-.-resultant 
coefficients 

MAPQ = element degree-of-freedom map 
NBT = NBM + NBB 

NQT = number of element intrinsic degrees of freedom 

X, X(T,J) = direction-I position relative to the 

element reference frame of node J 

NI 1 NI(n) = 1 , NJ(n) - j, the n*"^ element boundary 
NJ J connects element nodes i and j . 

TGEN directs the computation of the boundary work matrix T (NBT, NQT). 

The procedure incorporated in TGEN is outlined below: 

Step A. (TGEN line numbers 13 and 14) 

Set up MAPU(3,5) via data statement. For 
MAPUdjJ) = K, the K— boundary motion 
function IS associated with the Jill boundary 
stress resultant for computing boundary 
work. MAPU(2,J) and MAPU(3,J) pertain to 
warped boundary calculations and are not 
currently used. 

(J = IP, pg 36 and K = lU, pg 47) 

Step B. (TGEN line numbers 18 through 20) 

Initialize the T matrix to zero 
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Steps C through G are repeated for each boundary. 


Step C. (TGEN line numbers Z4 through 37) 

Set up the boundary geometry as; 

T1 = X location of boundary ongxn 
T3 = y location of boundary origin 
ZL = boundary length 
CS = cos G, 0 = angle between x and x 
SN = sin 0 


Steps D through G are repeated for each boundary stress resultant. 


Stress Resultant Identxfication; 


IF 

= 

1, 

Nu 

IP 


2, 

Nv 

IP 

= 

3, 

Mu 

IP 

= 


Mv 

IP 


5, 

^w 


membrane resultants 
active if NBM > 0. 


bendxng resultants 
active if NBB > 0. 


Step D. (TGEN line numbers 40 through 55) 

For membrane stress resultants (IP=1 and IP=2) 
TTGEN is called with a final argument of 1. For — 
bending resultants TTGEN is called with a final 
argument of 2. 

Calls to TTGEN, RGEN, and TCB result in the 
computation of arrays representing and T^^ 

of Eqs. 40 and 43. NMC represents X and C 

contains T , . ^ 

cb 


Step E. (TGEN line numbers 56 through 68) 

The boundary deformation function associated 
with stress resultant IP is set up by AGEN 
with the final argument lU. lU = MAPU(1,IP). 
NMA represents X^ (Eq. 41). The transformation 
matrix Tgq (Eq. 43) is stored in a mapped 
fashion according to arrays A, NQA, NQl, and 
IDQ. 
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Step F, (TGEN line number 70) 

Calling NCALNA results in the calculation 
of 


X 

( a 


ZL 

/ X X*^ dx (Eq. 42). 

J ^ " 


0 


X IS stored in the array XNCNA. 

C3. 


Step G. (TGEN line number 77) 

CXA performs the transformation 

t\ X T (Eq 43), 

cb ca aq ^ 

The result is added to the e^cisting 
T matrix. 
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TTGEN (PARENT ROUTINE = TGEN) 


The following inputs are supplied through the TTGEN calling sequence: 


IKIND = I, membrane stress resultants 

= 2, plate-bending stress^ resultants 


Tl 

T2 I X = Tl + T2 X 
T3 [ y = T3 + T4 X 
T4 


boundary line equation 
x-y = element ref, frame 
X = boundary ref- frame 


( 45 ) 


The following arrays are returned through the TTGEN calling sequence: 

TT = transformation matrix set up by 

TT6X3 or TT10X3 

NMC(1,NC) = representation of X (Eq. 40) 

NMC (2,NC) is not used. 

The procedure incorporated in TTGEN is outlined below. 

Step A. (TTGEN line numbers 4 through 6) 

For membrane stress resultants, TT6X3 
is called to generate TT(NX,NC). 

Step B, (TTGEN line numbers 8 through 10) 

For bending stress resultants, TT10X3 
is called to generate TT(NX,NC). 

Step C. (TTGEN line numbers 12 through 16) 

Set up NMC so that NMC (1,1) = n 

indicates that X (I) = x”^. 

c 
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IT6X3 (PARENT ROUTINE = TTGEN) 


The followxng inputs are supplied through the TT6X3 calling sequence. 


T1 

T2 

T3 

T4 


► 


boundary line equation, ?ee TTGEN 


The transformation matrix T is returned through the TT6X3 calling 
sequence. Where x and y refer to element reference frame axes and 
X refers to the boundary reference frame axis. 



xy 


T(6,3) 



or X = T X 
m c 


(46) 


The X^ vector corresponds fo the assumed membrane stress-resultant 
distribution presented in the K2D discussion 


K2D-39 



TT10X3 (PAEENT ROUTINE = TTGEN) 

The following inputs are supplied through the TT10X3 calling sequence: 


T1 

T2 

T3 

T4 


boundary line equatxons, see TTGEN 


The transformation matrix TT is returned through the TT10X3 calling 
sequence. Where x and y refer to element reference frame axes and x 
refers to the boundary reference frame axis, 



The X^ vector corresponds to the assumed bending stress-resultant 
distribution presented in the K2D discussion. 


K2D-40 



RGEN (PARENT ROUTINE = TGEN) 


The following inputs are supplied through the RGEN calling sequence* 


IP - 1. 

= 2, Nv 
- 3, Mu ^ 
= 4, Mv 
“ 5 , 

SN = sin 0 
CS = cos 0 . 


IP identifies the active 
boundary stress resultant. 


Define the orientation of the boundary 
relative to the element ref. frame. 


The array R is returned through the RGEN calling sequence. R is used 
to perform the transformation from element stress resultants to 
boundary stress resultants. The figures below indicate the transformation 
matrices that are required. Directions x and y refer to element reference 
frame directions. 



(boundary ref. frame axis) 


N 1 

u 

N 

V 


-SN^ -CS^ 2SN-CS 
-SN CS SN-CS (CS^-SN^) 


N 

y 

N 


u ^y. 


(48) 



X (boundary ref. frame axis) 



M 

X 

r M 1 
11 

= 

*SN CS 

-SN CS 

2 ? 
(CS -SN^) 

- 

m' 

xy 


= 

_-SN^ 

-CS^ 

-2CS SN _ 


M 


M 


M 


( 49 " 


= [-SN CS ] 


L^vj 


xyj 
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The R array is computed such that for 


IP=1 : Ny = R(3) 


IP=2: 


N^=R(3) Tn, 


IP=3 ; = R(5) 


, IP=4: 


= R(5) 


Mjj , and 
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TCB (PARENT ROUTINE = IGEN) 


Calculation of the transformation matrix (Eq. 43) is based on the 
following derivation. 


\ = (Cl C2 




1 


X 

x2 


a boundary stress (55) 

resultant in boundary 
ref. frame coordinates. 

See Eq 40. 


'si ■ 


®11 Ri2 


’ 1 ' 

9 

S2 

= 



X 


* 



■ 

y 


L^niJ 


-®NI,1 \ *** 


• 

• 



the assumed element (56) 

stress-resultant distribution 
(either membrane or bending) 
in element ref. frame 
coordinates b^ = 3— stress 
coefficient. See K2D. 


\ - (Ri R2 




51 

52 



transformation from element (57) 
stress resultants to 
boundary stress resultants. 

See RGEN. 


1 

X 

|y 


Zll Zl2 

Z 21 


• Z 


1,NC 


X 

x2l 




"NJ,NC-J 


transformation from (58) 

boundary coordinates to 
element coordinates. 

See TTGEN. 


Substituting Eqs. 58 and 56 into Eq. 57 yields the following 


(Rl R 2 ••• Rj^j) 


Rll ^1 R12 1^1 " ' 


“Zll Zi2--- 


‘ 1 ■ 

• 


Z,2i 


X 

» 


* 


x^ 

• 

B b ' ' " 




• 



7 7 




-NJ.l NJ,NC-‘ 



(59) 


- (Cj C2 ••• Cj^^) 


x2 
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( 60 ) 


Comparing like terms of the above equation (i.e like powers of x) yields 
the desired transformation matrix as: 


”ci ‘ 
C2 

= tcb(nc,nb) 

*bi ■ 

b2 

-Sc- 


- Sb- 


The following inputs are supplied through the TCB calling sequence. 


NI 

NJ 

B(NI,NJ) 

IB(NI,NJ) 


matrices characterizing the 
assumed stress-resultant 
distribution. These arrays 
originate in K2D. 


R(NI) = element stress resultant to boundary 
stress resultant transformation 
array. R is created in RGEN. 


NC 

Z(NJ,NC) 


transformation from boundary 
coordinates to element coordinates . 
Z is created in TTGEN. 


NB = ninnber of stress resultant coefficients 

used in the formulation 


The matrix TCB(NC,NB) is computed and returned through the TCB calling 
sequence 
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AGSN (PARENT ROUTINE = TGEN) 


Displacements and rotations relative to the boundary reference frame 
and relative to the element reference frame of boundary i-j end points 
are shown below. 


y 





For boundary x-j a deformation function, (Eq. 41), of the following 
form IS assumed: 


Uj^ = (1 5 ... ) 


ai 

a-2 


L^NA 


= polynomial wxth NA forms 


( 61 ) 


= a 
a 


The vector of polynomxal coefficients is expressed' in terms of the 
motions of nodes x and’ so that 


ai 

^2 

• 

• 

• 


1 

AB(1',NA,NQ1B) ' AB(2,NA,NQ1B 
1 

a «« 

t 

-i 

‘lit 

^m 

NQIB = 2 
xn thxs 
equation 

/na. 
















(62) 


Transformation from the element reference frame to the boundary 
reference frame yields 



D(NQ1B,NQ1) 


rnO- 

<ln 

NQl = 3 ' 


in thxs 

% 

equation 

1 



h = i or j . 


(63) 


For the above equatxon, the AGEN array IDQ would contain: 
IDQ(l) = n, IDQ(2) = o, and IDQ(3) = p. 
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Performing the transformation indicated by the previous equation 


yields 



i 

I 


L^Aj 




A(NA,NQA) 



‘n 



, where 



AB(1,NA,NQ1B) • D | AB(2,NA,KQ1B) 


“ 

D 


( 64 ) 


The following inputs are supplied through the AGEN calling sequence 


ZL = boundary length 


SN = sxn0| orientation of the boundary 

CS = cos 9 1 relative to the element ref. frame 


lU 


1, displacement in x direction 

2, displacement in y direction 

3, displacement in z direction 

4, rotation about x-axis 

5, rotation about y-axis 

6, rotation about z-axis 

7, for warped boundary (not currently used) 


The following arrays are returned through the AGEN calling sequence: 


A(NA,NQA) 

IDQ(NQl) 

NA 

NQA 

NQl 




> 


the meanings of these arrays are 
as indicated in the preceding 
discussion 


NMA(1,NA) 


representation of 

NMA(2,NA) IS not currently used 
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The procedure incorporated xn AGEN is outlxned below 


Step A. (AGEN line numbers 4 through 106) 

- ^ ^ ^ {1 

Set up NA, NQl, NQIB, IDQ, AB, and D 
according to'the xnput value of lU. 

lU - 1 (AGEN Ixne numbers 11 through 20) 

» T 

= displacement xn x-dxrection 

= aj + ag X 


’-1/L i 1/L 


f L^l 


[CS Sn] Jq^ 


lU = 2 (AGEN Ixne numbers 22 through 42) 

= displacement in y-direction 
-aj + a2X + a3 x^+ a^ 


-3/l2 -2/L I 3/L^ 
. 2/l3 1/l2 I -2/l3 
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lU - 3 (AGEN line numbers 44 through 53) 

Uj^ = displacement in z-dxrection 
= aj + a2 X 


ai 


1 

0 


"-X ' 
^3 



-1/L 

1/L 


M- 

^3 

= 

[SN 

-cs] 

q? 





h 

L^aJ 


lU = 4 (AGEN Ixne numbers 55 through 64) 
IL = rotation about the x-axis 




1 1 
1 

0 


% 



-1/L 1 

1/L 


M- 


?4 - [CS SN] 



lU = 5 


U, 


(AGEN line numbers 66 through 74) 
rotatxon about the y-axxs 


ai 



[l/L 

-1/l] 

q^ 

3 



qJ 

L 3-* 

[SN 

-CS ] 

r^r 





( 69 ) 


(70) 


(71) 


(72) 


(73) 


( 74 ) 
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lU = 6 (AGEN line numbers 76 through 95) 


U, = rotation about the z-axis 
k 


ai 

+ a2 

X + 33 x^ 



ai 


’ 0 1 1 

0 

» 1 

a2 


-6/L^ -4/L 1 

6/h^ -2/L 



6/l3 3 / 1 ? 1 

-6/l3 3/L^ 

rrh"" 

^2 


10 0 




5 ^. 


0 SN -CS 







„h 








^6 


4 


( 75 ) 


(76) 


III = 7 (AGEN line numbers 98 through 104) 

= displacement in the y-direction for 
warped membranes = a + a x 


ai 

' 1 1 
= 1 

° 1 

(77) 


_ -1/L 1 

1/^ J bl 


^2 



(78) 


Step B. (AGEN line numbers 107 through 111) 

Set up NMA(1,I), I = 1,NA. NMA 
characterizes the vector X , so that 
for ® 

X (J) = 

ci 

NMA(1,J) = k. 

NMA(1,1) = 0, 1. 

NMA(1,2) = 1, X 

NMA(1,NA) *= NA-1, 

NMA(2,NA) is not currently used. 
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Step C. (AGEN line numbers 112 and 113) 

Call ATD to perform the AB-to-A 
transformation. 
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ATP (PARENT ROUTINE » AGEN) 


The following inputs are supplied through the ATD callxng sequence- 


AB(2, NA, NQIB) 
D(NQ1B, NQl) 

NA 

NQIB 

NQl 




The contents of these arrays 
is discussed under AGEN. 


NQA = 2 * NQl 


ATD computes and returns the array A through the callxng sequence. A 
is confuted as 


A(NA,NQA) = J^AB(1,NA,NQ1B) * D 1 AB(2,NA,NQ1B) • pj . (79) 
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NCALNA (PARENT ROUTINE = TGEN) 


The following inputs are supplied through the NCALNA calling sequence 


XL 

= boundary length 

NMC, 

for X^(J) = x^, NMC(1,J) = k 


NMC(2,J) is not used 

NMA, 

for X^(J) = NMA(1,J) = k 


NMA(2,J) is not used 

NC 

= order of X 

c 

NA 

= order of X 

a 

NALFA' 

these arrays pertain to warped 

ALFA 

element boundaries and are not 

SIGN 

currently used 


NCALNA computes and returns through the calling sequence the array Z. 
Z is computed in accordance with Eq. 42 as 

XL 

Z(NC,NA) = f X X* dx 
yea 

0 

The procedure incorporated in NCALNA is outlined below 


Step A. (NCALNA line numbers 8 through 13) 


r -n - fXLr' ^ 

Z(I,J) = J X dx = , where 

N = NMC(1,I) + NMA(1,J) 

Since NCALNA is called several times for each 
element boundary, the integral of x^ may be 
required repeatedly. To avoid duplicate 
calculation, the common block /CK2D5/ containing 
IQX(IO) and QX(IO) is maintained, so that 

XL 

QX(N+1) = J dx. 

0 



If IQX(I) = 0, QX(I) has not been computed. 
If IQX(I) 0, QX(I) has been computed and 
is present in CK2D5. 


Step b. (NCALNA line numbers 15 through 24) 

This portion of the code pertains to warped 
boundary considerations and is not currently 
used. 
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CXA. (PARENT ROUTINE = TGEN) 

The following inputs are supplied through the CXA calling sequence: 


C(NC,NBC) = of Eq. 43, from TCB 

X(NC,NA) = of Eq. 43, from NCALNA 

A(NA,NQA) = of Eq. 43, from AGEN 


NA 1 
NQA 1 
NG [ 
NBC J 


dimensions of the above matrices 


NQT = number of element intrinsic degrees 
of freedom 

NBT = number of assumed stress resultant 
coefficients in the formulation 


INCB = row increment 
MAPQ(4,5) = degree-of-freedom map 

IIXKNQl) 1 

NQl 1 

NODES = 2 for current application 

NIJ(l) =11 this boundary connects nodes 
NIJ(2) = j J i and j 


The following arrays are returned through t^e CXA calling sequence: 

(83) 

(84) 

The procedure incorporated in CXA is outlined below: 

Step A. (CXA line numbers 9 through 13) 

Compute TMP(NBCjNA) = C X 


T (NBT, NQT) = T (NBT, NQT) + C X A 

A 

TMP(NBC,NA)= C X = temporary array 
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Step B. (CXA l^.ne numbers 15 through 27) 


Compute T(NBT,NQT) = TMP(NBC,NA) * A 

The row location in T of the product IMP -A 
is controlled by the increment INCB. The 
column location is controlled by the three 
maps IDQ, and NIJ. 


NIJ(l) NIJ(ND) = NODE NIJ(NOPES) 

- - ■ ^ ^ 


1 

*■ 1 

1 ' 1 

1 

i 

|.r2 '*- IQ •• NQ'i 


j nr] 

■ 

2 









A(NA,NQA) ’ 


ND = 1 

« * » 



« * « 

ND = NODES 


» 

NA 


- 



- 


_ 



} 

IDQ (IQ) = IIQ 


MAPQ(NODE,IIQ) = JT 

The product of IMP with this 
column of A is added to 
column JT of T. 


Note that NODES = 2 fqr the current 
application. 
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SSTM (PARENT ROUTINE = K2D) 


The following inputs are supplied through the SSTM calling sequence: 

H(NHD,NHD) = Inverse of the element 

strain energy matrix 

T(NHD,NQ) = boundary work matrix 

dimensions of the above matrices 
NQ 

NB = NHD in the current application 

The following arrays are returned through the SSTM calling sequence: 

STM(NB,NQ) = element stress matrix 
S(l) = element stiffness matrix 


The following procedure is incorporated in SSTM: 
Step A , (SSTM line numbers 3 through 7) 
Compute STM = HT 


Step B. (SSTM line numbers 8 through 14) 

^ it 

Compute S = T HT = T STM. S is 
stored by column such that 


S = 


Sj S2 Sif 

S3 S5 
Se 

« 


1 H = NQ AtL ) 



( 86 ) 


(87) 
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PROGRAM LISTING 


The routines comprising the K2D element formulation are listed, in 
the order shown in the table of contents on this page and the following 
18 pages. 


K2D 


SUBROUTINE K2D ( XX, S^STM> NS ICES, NBM^NBB# ISHELL, NOT, X 34, 1X3^, I ERR ) 

DIMENSION A(1),X(2#4),XX{2,1) 

♦, IBH{3,6),BM{3,6) ,XYM(6),IXYM(2,6>) 

♦ , IBB(5,10),BB(5,10),XY8( 10),IXYB(2, 10) 

MAPH(4,5),MAPB(4, 5) , MA PQ ( 4, 5 ) , NQH ( 2 ), NQ6 { 2 ) , NI ( 4 ) , N J I 4 ) 

W(IOOO) 

COMMON/ WARP/ I WAR P, NAL F A ( 4 ) , ALF A ( 5, 4 ) , NGAHA ( 4 ) , GAM A ( b , 4 ) 


9 


C0MM0N/£RRK2D/AERR,NERR 



10 


C0MM0N/CK2D1G/JUMP(1C) 



11 

C + 

LABELED COMMON BLOCKS USED WITHIN 

K2D ROUTINES ARE 

12 

C* 

/CK2D1/DAT(11) 



13 

C + 

/CK2D2/DAT(105) 



14 

C* 

/CK2D3/DAT(3) 



lb 

C + 

/CK204/DAT(49) 



16 

c* 

/CK2Db/DAT(20) 



17 

C + 

/CK2010/DAT( 10) 



16 

c 




19 


DATA ( (ISMU, J), J»1,6),I«1, 3), 

( (BMC 1,J) 

,J*1,6), I 

20 


* 3^ 6, 4^ 10, 6, 0, 



2i 


♦ 2, 5, 7, 9, 10, 0, 



22 


* 1, 7, 6, 0, 0, 10, 



23 


+ 1«, 1», 1«, 1«, 1», o«. 



24 


^ 1», 1«, 1», 1», 1«, 0*^ 


* 

2b 


lo,*"l«,~l», 0», 0»,”2» / 



26 

c 




27 


DATA (aXYM(I,J), J«1,6),I»1,2) 

, (XYH(I), 

I«l,6) / 

26 


^ 0, 1, 0, 2, 0, 1, 



29 


♦ 0, 0, 1, 0, 2, 1, 



30 


* 1», 1», !•, 1«, 1*, !• / 



31 

c 




32 


DATA ((IBB(l,J),J«l,10),I»i,5),((tJB(I,J),U»l,10) 

33 


♦ 1, 0, 6, 0, 4, 0, 14, 

12, 10, 

0, 

34 


* 2, 0, b, 0, 7, 0, 13, 

lb, 11, 

0, 

3 b 


■S' 3, 0, 9, 0, 8, 0, 0, 

0, 14, 15, 

36 


♦ 6, 6, 14, 15, 10, 0, 0, 

0, 0, 

0, 

37 


^ 7, 9, 11, 0, lb, 14, 0, 

0, 0, 

0, 

36 


^ lo, 0*, 1«, Oe, 1», 0*, 1*, 

1 e , 1 « , 

Op, 

39 


♦ Ic, 0«, 1», 0», i», 0«, 1«, 

1 • , 1 •, 

Op, 

40 


1«, 0*, I. • , 0», 1», 0«, 0», 

0 • , 1 • , 

1 ® , 

41 


♦ lo,”!*, I*,”!,, 1», 0», 0», 

0 i , 0 , , 

Op, 

42 


♦ I*,"”!*, 1», 0«, I*,"”!*, 0«, 

0® , 0 », 

0. / 

43 

c 




44 


DATA ((1XYB(I,J),J»1,10),I»1,2),(XYB(I) 

,I»1,10) 

4b 


♦ 0, 0, 1, 1, 0, 0, 2, 

0, 1, 

1, 

46 


♦ 0, 0, 0, 0, 1, 1, 0, 

2, 1, 

1, 
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47 


♦ 

1 o> 1 • 9 

laf 

1 o ^ 

Itf laf la f 

laf 

48 

C 








49 


DATA 

(N0M( I) 

= 1,2); (NQBdJiI* 

XfZ) 

50 

c 








5i 


DATA 

( (MAPM(Ij) J J 

•1/5)#I*1>4) 

/ 

52 


* 

0, 

0, 

Of 

0, 

Of 


53 


* 

If 

0, 

Of 

0^ 

Of 


54 



If 

3f 

Of 

0^ 

Of 


55 



4# 

5» 

Of 

Of 

0 / 


56 



DATA 

( (HAPBd# 


*l,5)d»l#4) 

/ 

57 



Of 

0^ 

Of 

Of 

Of 


5b 



Of 

0, 

If 

Zf 

3f 


59 


❖ 

Of 

Of 

4# 

5f 

bf 


60 



Op 

Of 

7f 

8f 

9 f 



KW“1000 

iWARP-0 

A£RR«4H 

NERR«0 

DO 600 I«l>2 
DO 600 J»l#3 
xa#j)“XX(i» j) 

JaNSIDES 

X(lp^)«XX(liNSiOES) 

X(2#^)=XX{2^NSI0ES) 

NB0Tri«O 

IF{NBM«GT.OeAND.NBB,GTeO) 

NBT*>NBM+NB8 


N80TH»1 


MA«1 


NSl°NSI0eS-2 

NMQoO 

NBQ«0 

IF{NBM,GToO) NHQ®NQM(NS1) 

IFCNBBoGTeO) NBQ«NQ8(NS1) 

NXM = 6 
MXB°10 
HA“0 
MB^O 

XFlKBH.GTeO) MA«1 
IF(N8BoGT*0> MB=1 
DO 100 INODE=l#NSIDcS 
DO 100 
HG=0 

IF(HAP?HINDDE> JJ.EQ.O) GOTO 500 
MQ»MA««APM{lNaDE# J) 

500 IFCMAPBUNODEW) .EQ.O) GOTO ^00 
MQ“MQ+MB«(rtAPB( INODEi J)tNHQ) 

^00 MAPQdNODE, J)»flQ 
iOO CONTINUE 

IFtNBSo&TeO) GOTO 200 
IF(IX34.£Q,0.0R.NSID£S,EQo3) GOTO 200 
NMQ«6 

MAPQ(^j.3)«6 
200 NQToNHQ+NBQ 

IHal 

NHoNBT«NBT 


GOTO 
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lOi! 

104 : 

105SC 

106: 

107: 

10b: 

109: 150 

110:C 
111 : 

112 : 

113: 

114: 

116: 

116: 

il7: 300 

llb:C 

119: 

120: 

121: 250 

122:C 

123: 

124:C 

125: 

126: 

127: 350 

128:C 

129: 

130; 1000 
13i: 

132: 2000 
133: 

134: 

135: 

136; 


KW»KW-NH 

IT-NH+l 

IF(4UMP(i).N£,0) GOTO 150 

CALL HGEN{NBM,NBB>W(1H)> (IT)# KWjX# A,XYH#IXyM>XYB»IXYB^ IBM, BM 
IBB#8B#NXM#NXB#ISHELL#NBT#I6RR) 

CONTINUE 

NT»NBT*NQT 

IF(NT*GT.KW) GOTO 1000 
Nia)«l 
NJ(1)»NSIDES 
DO 300 I»2#NSIUES 
NI (I )»NJ( I-l) 

NJ(I)«NI(I)-1 

IF{ JUMP(2),NE.0) GOTO 250 

CALL TGENdBM, BM# NBM# IBB# BB# N6B# NS IDES# MAPCi, W ( IT ) # NB T# NOT# X# M# NO ) 
CONTINUE 

IF(1X34.6T.0) CALL WARPT ( W ( I T ) # HAPQ# NfaM # 0# NBT #NQT# X# X 34 ) 

IF(OUMP(4).NE.O) GOTO 350 

CALL SSTM(W(IH)#W(1T)#S#STM#NBT#NBT#NQT) 

CONTINUE 

RETURN 
N£EO=NT«Kw 
WRIT£(6#2000) NEED 

FORMATC SlHO’i^*^ ADDITIONAL CORE REQUIRED FOR ELEMENT K# NEED ■ 

AERR-4HC0RE 

NERR-NEED 

RETURN 

END 


HGEN 


i: 

2 : 

3: 

4: 

5: 

6 : 

7; 

b: 

9: 

lO: 

li: 

12 : 

13: 

14;C 


SUBROUTINE HGEN ( N8M# NBB# H# W# KW# X# AA# XYM# IXYM# XYB# I XYB# I BM# BM# I BB 
♦# BB#NXM#NXB#ISHELL#NBT#IERR> 

DIMENSION X(2#4)#W(1)#A(3#3#3)#H(NBT#NBT)#AA(1) 
C0MN0N/ERRK2U/AERR#NERR 
C0MM0N/CK2D4/IQT(7#7)#QT(7# 7) 

CDMM0N/CK2D1/XX (4)#YY(4)#NXD#NYD#DUH 
COMMON/CK2D10/ JUMP(iO) 

DATA NC#NSM#NSB / 3#3#5 / 

DO 160 I»l#49 
160 IQT(I#l)-0 
DO 260 I*l#4 
XX(I)»X(1#I) 

260 YY(I)»X(2#I> 
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lb: 
16J 
17: 
lasC 
19: 
20 : 
21 : 
22 : 
23: 
2^ : 
25: 
26: 
27: 
2 &: 
29: 
30: 
31: 
32: 
33: 
34: 
35; 
3b : 
37: 
38: 
39 : 
40; 
4i: 
42 : 
43: 
44: 
45: 
46: 
47: 
46: 
49: 
50;C 
51: 
52: 
53: 
54; 
55: 
56 ; 
57; 
58: 
59: 
60 : 
61 : 
62: 
63: 
64: 
65: 
66 : 
67: 
68 : 
69: 


DO 800 I°1^NBT 
00 800 J«1^NBT 
800 H(IiJ)nOoO 

NBOTHaO 

IF(NBM,GT*0.AND,NBBoGT,0) NBDTH-1 
00 600 I»1^27 
600 A(I^lj)l)sOoO 

lF(NBMcL6.0) GOTO 400 
IF(NBH.GT»1) GOTO 500 
A(3>-3,l)»loO/AA(l) 

GOTO 400 
500 K»0 

DO 700 J«li>3 
DO 700 
K»K+1 

A(I> Jpl)»AA(K) 

700 A( I> l)oAA(K) 

400 IF{NBBcLEoO) GOTO 350 
K°6 

IF(NBMoLE<.0) K®0 
00 900 J«1 j 3 
IF(NBMo6T«0) K«K+3 
DO 900 I°1#J 
KbK+1 

A(Ij>»)# 3)oAA(K) 

900 A( J> Ip3)»AA(K) 

350 IF (NBOTH.LE.O) GOTO 150 
K°6 

DO 250 J«l>3 
KoK«J-l 
DO 250 Iol>3 
K=K-H 

250 Adi. Ji.2) = AA(K> 

150 CONTINUE 

NEXT'^1 

IFtNBHoLEoO) GOTO 100 

iPXf1°NEXT 

IPYN°IPXH+NC^NBH 

IPM°IPYtl4-NC«N8M 

NEXT“IPh+NC*NBM 

IFCNEXToGToKfe) GOTO 1000 

CALL HGENDCNBN^ 0>NBT#NC#A(1>1/1)>XYM^ IXYM>W<IPXM) ,W(IPYM) 

*t W(IPM)#H^NXN^ IBH>BM,NSM#0) 

100 IF(N8B,LEoO) GOTO 300 
IPXB»NEXT 
IPVB»1PXB4-NC*NB8 
IPB®IPYB+NC*NBB 
NEXT“IPB+NC*NBB 
IF(N£XT«GT*K!s) GOTO 1000 

CALL HGEN0(NBBi NBM>NBT»NC^ A(l#l> 3)# XYB, IXYB^W(IPXB)^ W (IPYB) 
*f W(IPB),H/NXB,IBB,BB,NSB^O) 

IF(ISHELL«EQ,0) GOTO 300 

CALL HNBGEN(NB MjNB6jNC»NCj A(l^ 1> 2 ) p W (1 PXN ) > W ( 1 PY H ) j kJ ( IPH ) 
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70S W(IPX6),W(iPYB)>W(IP8)#H>NBT) 

7i: 300 CONTINUE 

72: IF( JUMP(3).NE.O) RETURN 

73: CALL SYMVRT(H>N8T,IS1NG#N8T) 

7A: IFdSING.EO.O) RETURN 

75s WRITE(6/2000> ISING 

76S 2000 FORHAT{ A3H0*** WARNING# H MATRIX IS SINGULAR# ISING » 14 ) 

77: IERR-2 

785 RETURN 

79s 1000 WRITF(6#3000) KW#NEXT 

80: 3000 FORHAK 38H0*«»«' INSUFFICIENT CORE FOR COMPUTING H / 

8is ♦ 110# 9H PROVIDED# liO# 9H REQUIRED ) 

82: IERR»1 

83: A£RR=4HC0RE 

84: NERR»NEXT 

85: RETURN 

86: end 


HGEND 


1: SUBROUTINE HGEND ( NB# INC B# NBT# NC# C# XY# IX Y# I PX > I PY# P# H# NX# I B# B #NS 

2: *, INCC) 

3s C0MH0N/CK2D4/IQT(7#7)#QT(7#7) 

4; C0MMQN/CK201/XX{4)# YY(4),MX#MY#Z 

5: DIMENSION C(NC#NC)#IPX(NC#NB)#IPYINC#NB)#?(NC#NB) 

65 +# H(NBT#NBT)#XY(NX)#IXY(2#NX) 

7: ♦# IB{NS#NX)#B(NS#NX) 

8:C 

9: 2ER0alo0E~10 

lose 

11 : 

12s 

13 s 
14S 

15; 400 

16: 

17s 
18 : 

19 : 

20: 

21 : 

22: 

23; 

24: 

25: 

26: 500 

27 • C 

2bs DO 100 I»1#NB 

29: IH°I+INCB 

305 DO 100 J»1#I 


DO 400 I=1#NC 

00 400 J«»1#NB 

IPXd# J)»0 

IPY(I# J)»0 

Pd# J)»0.0 

DO 500 IC»1pNC 

IIC°IC+INCC 

DO 500 IX»1#NX 

IIBoIB(HC#IX) 

IFdlB.LE.O) GOTO 500 

IFdIB.GT.NB) GOTO 500 

BIJ°B{IIC#IX) 

P(IC#IIB)»BIJ+XYdX) 

IPX{IC#II8)*1XY(1#IX) 

IPYdC#IIB) = IXY(2# IX) 

CONTINUE 
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3i: jH»J+INCB 

32: DO 200 IC»1>NC 

33: DO 200 JC«1,NC 

34: PC»P<IC#I)*C(1C#0C)*P{JC# J) 

3i)» If (ABS ( PC) .LToZERO) GOTO 200 

36: MX»IPX(ICj>I)+IPX{ JC#J) 

37: MY-IPY(IC,I)+IPY(JCj J) 

3tJ: HMX-MX + 1 

39; rtMY»hY+l 

40; IF(IQT(HMX^M«Y) .EO.l) GOTO 600 

41: CALL_ITQUAD 

42; I0T(M«XpMMY)»1 

43; QT(MMX>MMY)«Z 

44; 600 H(IH, JH)+PC*OT(MMX,HMY) 

4b; 200 CONTINUE 

46: H( JH,IH)«H(Ih# JH) 

47; 100 CONTINUE 

4b;C 

49: RETURN 

bO: END 


HMBGEN 


i; SUBROUTINE HhBGhN(NBM,NBB^NCM>NCB,C,IPXH, IPYM>PM>IPXB,IPYB, PB 

2 ; *, H#NBT ) 

3: OINENSION C ( NCM, NC B ) i IPXM ( NC M, N8M ) # I P YN ( NC NBM ) 

4; PH(NCH,NeN)#lPXB(NCB,NBB), IPYB(NCB,NBB)# PB(NCB,NBB) 

5; HJNBTpNBT) 

6: C0MH0NyCK2D4/IQT{7>7)j>QT{ 7#7) 

7: C0MM0N/CK2D1/XX{4J,YY{4),NX>MY»Z 

B:C 

9: ZfcRD«'1.0E“10 

io;c 

11 :C 

12: DO 100 1«1,NBH 

13; IH«1 

14; DO 100 0=1^NB8 

15; JH-J+NBM 

16; DO 200 IC«1,NCM 

17: DO 200 JC»1>NCB 

18 ; PC«PM( IC, 1 )«C( IC/ JC)*PB ( JC^ J) 

19: IF(ABS( PC ) oLT.ZERC) GOTO 200 

20; NX*IPXM(ICj I) + IPXB( JCW) 

2i: MY»IPYM{IC,I)+IPYB( JC^ J) 

22: NNX-NX+1 

23: flMY»«y+l 

24; IF(IQT(NNX#hriY).EQ.l) GOTO 600 

25; CAU.JTQUAD 
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26; 


IQT(NNX#MMY)“1 

27: 


UT{ NNXj>«MY)-Z 

28; 

600 

ri{ IH;. JH ) + PC«QT(NNX»«rtY) 

29: 

200 

CONTINUE 

30: 

100 

H{ 

31: 


RETURN 

32: 


END 




SYMVRT 

i 


SUBROUTINE SYKVRT ( A, N, IS ING, NHD ) 

2 

C 

SYHETRIC INVERSION SUBROUTINE 

3 


DIMENSION A(NHD,NHD)>6{30)>T(30),ERR (2) 

h 


CGMHDN/ERRK2D/AERR>NERR 

5 


DATA ERR(l)^ ERR(2)/4HISNG>4HNEGD/ 

6 


IF(NoGTol) GOTO 100 

7 


A(li»l)»l,0/A(l#l) 

6 


RETURN 

V 

100 

CONTINUE 

10 


IX“N 

11 


ISIN6-0 

12 


DO 18 I0*1*IX 

18 


T(ID)»1. 

14 


IDM1=ID-1 

15 


IDP1«ID+1 

16 


IF(IDH1)10^14,10 

17 

10 

DO 11 J«IO#IX 

18 


00 11 K*l^IDMi 

19 

11 

A(IDW)«A(ID, J)-A(K>ID)*A(K, JM'TIK) 

20: 

14 

IF{A(ID#I0))13yl2>15 

21 

12 

ISING^ID 

22 


NERR-ID 

23 


AERR«ERR( 1) 

24 


GO TO 30 

25 

13 

T{ID)s-l, 

26 


NERR«NERR+1 

27 


AERR»ERR(2) 

28 


AdD^ID) « ABS ( A(ID> ID) ) 

29 

15 

A(ID»ID)«SQRT <A(ID#ID))*T{ID) 

30 


IF(IDPl-IX) 16, 16,19 

31 

16 

00 17 J»I0P1,IX 

32 

17 

ACID, J )«A(IO, J)/A(ID^ ID) 

33 


A(ID,ID)«A(ID,ID)*T(ID) 

34 

18 

CONTINUE 

35 

19 

A(IX,IX)bA(IX,IX)*T( IX) 

36 


1XM1=IX-1 

37 


DO 24 l»l,IXhl 

38 


I1°I+1 
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39 

40 

41 
4 ^ 

43 

44 

45 

46 

47 
46 

49 

50 

51 

52 

53 

54 

55 

56 

57 
56 

59 

60 
61 
62 

63 

64 

65 

66 
67 


00 20 J«I1^IX 

20 A{I,4)«A(I^ J)/A(I,I) 

24 A{I#I)-T{I)/A{I^I> 

A(IX;.1X)»T{IX)/A(IX#IX) 

DO 23 I«2»IX 
K»I+1 

00 23 

IF(K-IX)21>23,23 

21 00 22 J-K^IX 

22 A(Li>J)“AlL; J)-A(I^ J)*A(L/I) 

23 A(L,I)=-A{L,i)*A(i,l) 

DO 27 I-IWX 

DO 29 K«I>IX 

29 B(I)«B(I)-s-AU,K)**2*T{K) 

M“I + 1 

IF(H-*IX)26>26,32 

26 DO 27 J“MjIX 
A{ 

DO 27 K«J#IX 

27 A(J^I)»A(J,I)+A(I»K)»A(J;K)*T(K) 
32 DO 28 I»1#IX 

28 A(I,I)«B(I) 

DO 35 J«1>N 

00 35 I«1^N 

35 A( 

30 RETURN 
END 


ITQUAD 


l: SUBROUTINE ITQUAD 

2: C0MM0N/CK201/ X( 4)# Y (4)^ N#H>AXNYM 

3: COHHON/CK202/ND»MD#MMDj CC(50),C(50),T,SUM 

4; C0MHQN/CK2D3/AR{3) 

5J C0MMDN/ZR0K2D/2ER0 

6: ND«N 

7: DO 100 I»1^3 

8: 100 AR{I)=0.0 

9:C* INTEGRATE ZONE 1 
10; X3hX4»X(3)-Xi4) 

ll; XL-X3MX4 

12: IF (X(2J .GT.XU XL»X(2) 

13: XM«1,0/XL 

14: IF(ABS(X3ttX4*XM) •LT.ZERO) GOTO 200 

15: A«<Y(4) 

16; B-(Y(3)-Y(4) )/X3MX4 

17: A«“1.0 

16: rtO-M+1 


K2D-65 



19: IF{ABS(8).GT.ZERD) GOTO 300 

20: MD“0 

2i: C(l)«loO 

22 : 

23: GOTO 700 

24:C* EXPAND THE Y BINOMIAL 

25: 300 AMP2°MD+2 

26: ADB^A/B 

27: C{1)»B«*MD 

28*. M1»MD+1 

29: DO 150 I«2^M1 

30: AI1«»I 

3i: AI-olAHPZ-AIl)/ (AIi-1.0) 

32: C( I)°Ca-l) + AI'*'ADB 

33: 150 CONTINUE 

3<f:C* EXPAND THE X BINOMIAL 
35: 700 AMP2-»N+2 

36: ADBoX(4) 

37: CC(1)»1.0 

38; IF(N.EQ.O) GOTO A-50 

39: Nl-N+1 

AO: DO 250 I«2#N1 

Ai: AI1°I 

A2: AI=(AMP2-AI1)/ (AI1-1»0) 

43: CC(I)«CC(I“1)«AI*ADB 

44: 250 CONTINUE 

45: 450 MHD»HD+1 

46: T»X3HX4 

47; CALL DSUM 

48: AR(1)«SUM 

49; AR (1)»AR(1)*AM 

50;C* INTEGRATE ZONE 2 

5i: 200 IF (ABS (X(4)+XH) .LT.ZEKO) GOTO 400 

52: ANM«N+M+2 

■53: AR(2)»{X(4)1'«(N+M+2) )«:'(Y{4)/X{ 4) >**(M+1 ) 

54: AR(2)°AR{ 2)/ANM 

55:C* INTEGRATE ZONE 3 
56: 400 X3MX2oX(3)-X( 2) 

57: IF(ABS{X3MX2>»XH> oLT»ZERO) GOTO 600 

58:C+ EXPAND THE X BINOMIAL 
59: AMP20N+2 

6O: ADB°X(2) 

6l: CCU)°loO 

62: iF<N»EQoO) GOTO 550 

63: N1°N-H 

64: DO 350 l«2j>Nl 

65; AI1«I 

66: AI»(AMP2-AI1 )/ (AI1“1*0) 

67: CCa)°CC(I-3 )+AI*AOB 

68: 350 CONTINUE 

69: 550 C{1)“1.0 

70: HD°M+1 

7i: MhO«l 

72: T“X3MX2 

73: CALL DSUM 
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74:C* sum 3 ZONES 


75: 


AR(3)»- SUM*{ Y(3)/X3MX2)**(M+1) 

76: 

600 

CONTINUE 

77: 


AXNYM«0 

78: 


DO 800 I«i^3 

79: 

800 

AXNYM«AXNYM+AR( I) 

80: 


AMl-H+1 

81 : 


AXNYH-AXNYM/AMI 

82: 


RETURN 

83: 


END 


DSUM 


1: 


SUBROUTINE DSUM 

2: 


C0MM0N/CK2D2/N,My MH>B(50)#C (50),T^SUM 

3: 


SUM«0. 

4: 


Nl-N+1 

6: 


M1«MM 

6 : 


IB»1 

7: 


NE«1 

8: 


NM«N+M+1 

9: 


lE-NM 

10: 


DO 100 ITERM-I^NM 

in 


NE-NE+1 

12 : 


Jl-NE-IB 

13: 


IF{J1.GT.M1) IB-IB+1 

141 


IF(IB.GT.Nl) GOTO 100 

15: 


lE-IE-1 

16 : 


D"0« 

17: 


DO 200 I»IB#M 

lb: 


J = NE--I 

19: 


IF(J.EO.O) GOTO 300 

20: 

200 

0»D+6(I)*C( J) 

21: 

300 

AE1«IE+1 

22: 


IEE-IE+1 

23: 


5UM-SUM+(D«T*+I££)/AE1 

24: 

100 

CONTINUE 

25: 


RETURN 

26; 


END 
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TGEN 


i: SUBROUTINE TGEN ( NBHj IB B # BB> NBB» NSI DES> MAPQ 

2: *, T#NBT>NQT,X,NI# NJ) 

3: DIMENSION IBM(3,1)^ BM(3,l)iIBB{5,l),BB(5il) 

4: ♦, T(N8T>NQT)>X{2>1)#NI(1)#NJ(1) 

bi *, NMA(2#9),NMC (2#9)>TT(100), XNCNA(9>9) 

6: ♦, R(5)^ A{9/8)>IDQ(4)^ HAPU{3>5),NIJ(2) 

7J ♦, C( 150),TEMP(250) 

B:C 

9; COMMON/WARP/ I W ARP^ NALFA ( 4 ) > AL FA ( 5^ 4 ) / NGAMA ( 4 ) , GAMA ( 5, 4 ) 

lOJ C0HH0N/CK2Di0/JUMP( 10) 

li: C0MM0N/CK2D5/IQX( 10)/QX(10J 

12:C 

13: data ( (MAPU(I,J), 1-1,3) W»1>5) / 3, 0,-2, 1, 2, 0, 6, 0,-5, 

14: ♦ iff 5, 0, 2,-1, 3 / 

15:C 

16: TH-1.0 

17:C 

la: DO 100 X«1,NBT 

19: DO 100 J-1,NQT 

■20: 100 T(I,J)-0.0 

2i: DO 200 IS I D£ -1, NS IDES 

22: DO 160 1-1,10 

23: 160 XQX{I)-0 

24: II-NUISIDE) 

25: JJ-NJ(ISIDE) 

26: T1-X{1,II) 

27: T3»X(2,II) 

2d: ZL-0,0 

29: 00 300 1-1,2 

30: 300 ZL-ZL+(X{I,JJ)-X(I,lI))++2 

31: ZL-SQRT(ZL) 

32: CS«(X{1,JJ)-X{1,II))/ZL 

33: SN-(X(2,JJ)-X(2,II) )/ZL 

34:C 

35: NOOES-2 

36: NIJ{1)-II 

37: NIJ(2)»JJ 

38 :C 

39: DO 500 IP-1,5 

40: IF( JUMP(5),NE.0) GOTO 700 

41: IF(IP,GT.2) GOTO 600 

42: IF(NBH.EQ.O) GOTO 500 

43: NBC-NBM 

44: IF(IP«E0.1) CALI TTGEN{TT,Ti,CS,T3,SN,NX,NC,N«C, 1) 

45: INCB-0 

46: CALL RGEN ( R, I P, SN, CS ) 

47: CALL TC 8 ( 3, NX, R, TT, 6M, IB M, C, NC, NBM ) 

46: GOTO 700 

49: 600 IF{NBB.EQ,0) GOTO 500 

50: NBC-NBB 
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51: IF(IP,EQ,3) CALL TT6EN ( TT> T 1, CS^ T3, S N# NX, NC >NMC > 2 ) 

52J INCB-NBM 

53: CALL R6EN (R, XP,SN,CS ) 

5<ts CALL TCB{ 5,NX,R,TT,BB,IBB,C,NC,NBB) 

55: 700 CONTINUE 

56* DO 800 IUU«1,3 

57: IU«HAPU(IUU, IP) 

58: IFdUU.EQol) GOTO 900 

59: IF(IWARP.EQ.O) GOTO 800 

60: IF(IUU.EQ«2.AND.NALFA(IS10E}.EQ.0} GOTO 800 

61: IF ( IUU.£Q.3.AND.NGAMA(ISIDE) .EQ.O ) GOTO 800 

62: IF(IU.EQ.O) GOTO 800 

63: iFCIUoGT.O) SIGN-+1.0 

64: IF(IU.LT.O) SIGN — 1.0 

65: IU»IA8S(IU) 

66: IF( IP.LT.3.AND.NBB.EQ.0) IU-7 

67: 900 IF( JUMP(6) .NE.O) GOTO 150 

68: CALL AGEN ( A, NNA, NA, NQA, NQl, IDQ, Z L, CS, SN, lU ) 

69: 150 IF(JUHP(7).NE.0) GOTO 250 

70; IF(IUU.EQol) CALL NC ALNA ( ZL, XNC NA, NMC, NMA, NC , NA, 

71: ♦ 0,ALFA,SIGN) 

72: IFUUU.EQ.2) CALL NCALNA( ZL,XNCNA,NMC,NHA,NC, NA, 

73: ♦ NALFAaSIDE),ALFA(l,ISIDE),SIGN) 

74: IFdUU.EQ.3) CALL NCALNA { ZL, XNC NA, nNC# NMA, NC, NA, 

75: ♦ NGAMA(ISIDe),GAMAd,ISIDE>,SIGN) 

76: 250 IF.{JUMP(8) .NE.O), GOTO 800 

77: CaLl CXA(C,XNCNA, a, T, TEMP, NA, Nq’a, NOT, NC , NBC, NBT, INC B, MAPQ 

7b: ♦, IDQ, NODES, NIJ,NQi) 

79: 800 CONTINUE 

60: 500 CONTINUE 

8i: 200 CONTINUE 

82: RETURN 

83: end 


TTGEN 


i: SUBROUTINE TTGEN ( TT,T1, T2, T3, T4, NX, NC, NMC, IKIND ) 

2: DIMENSION TT(l),NMCt2,l) 

3: GOTOdO, 20), IKIND 

10 NX-6 

5: NC»3 

6: CALL TT6X3 ( T 1, T2,T3, T4, TT ) 

7: GOTO 100 

8: 20 NX-10 

9: NC-3 

10: CALL TT10X3(T1,T2,T3,T4,TT) 

11: 100 CONTINUE 

12: NXP — 1 

13: DU 30 1-1, NC 

14: NXP-NXP+1 
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15: NMC<1>I)-NXP 

16: 30 NMC(2,I)«0 

17: RETURN 

Itt: END 


TT6X-3 


i: SUBROUTINE TT6X3 (T T2, T3, T4, T ) 

2: DIMENSION T(6#3) 

i: 00 150 1-1,6 

4; DO 150 J-1,3 

5: 150 T(I,J)»0.0 

6: T(l,l)-1«0 

7: T(2>1)-T1 

b: T(3,1)»T3 

9: T(4,1)-T1**2 

lO: T(5,1)«T3**2 

li: T(6,1)«T1*T3 

12: T(2,2)»T2 

13: T(3,2)«T4 

14: T(4,2)-2.0*T1+T2 

15: T(5,2)-2.0>:'T3+T4 

16: T(6,2)-T1*T4+T2*T3 

17: T(4,3)«T2+*2 

lb: T(5,3)-T4**2 

19: T{6,3)«T2+T4 

20: RETURN 

21 : END 


1710X3 


i: SUBROUTINE TT10X3 ( Tl, T2, T 3, T4, TT ) 

2: DIMENSION TTC10,3) 

3: DO 200 1-1,10 

4: 00 200 J-1,3 

5: 200 TT(I,J)-0.0 

6: TT{1,1)-1.0 

7: TT(3,1J-T1 

8: TT(3,2)»T2 

9: TT(5,1)-T3 

lO: TT(5,2)-T4 

ll: TT(7,1)-Ti**2 

12: TT(7,2)-2.*T1*T2 

13: TT(7,3)«T2**2 

14: TT(8,1)-T3**2 
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Lbt 

16: 

17: 

18: 

19: 

20 : 

2i: 

22 : 

23: 

24 : 

25: 

26; 


TT(8#2)»2.*T3*T4 

TT(8,3) «T4+*2 

TT(9,1)»T1+T3 

TT(9,2)«T1*T4+T2*T3 

TT(9,3)«T2*T4 

00 300 1«1^4 

II»2*I 

1F(I.EQ.4) II-IO 
DO 300 J-l#3 
300 TT(1I, J)-TT(II-1, J> 
RETURN 
END 


RqEN 


1 

2 

3 

4 

5 

6 

7 

8 
9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 

27 

28 

29 

30 


SUBROUTINE RGE N { R, I P, SN, C S > 
DIMENSION R(l) 
G0T0(10^20,30>40>50)^ IP 
10 R(1)«-SN**2 
R(2)«“CS*+2 
R (3 )»2.0*SN*CS 
GOTO 100 
20 RU)»-SN*CS 
R(2)-“R (1) 

R (3)««CS**2-SN**2 
GOTO 100 
30 R{1)»SN*CS 
R(2)*“SN*CS 
R(3)»CS**2-SN**2 
R (4)«0.0 
R(5)»0.0 
GOTO 100 
40 R(l)— SN+ + 2 
R{2)»- CS**2 
R (3)— 2. + SN*CS 
R(4)»0.0 
R(5)-0.0 
GOTO 100 
50 DO 200 I«1^3 
200 R{I)»0.0 
R (4)—SN 
R(5)»CS 
100 CONTINUE 
RETURN 
END 
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TCB 


l: 

SUBROUTINE TC B < NI, N R<# Z# 6, IB# XCB,NC, NB ) 

2: 

DIMENSION R( 1)#Z(N0#NC)#B<NI#NJ)#IB( NI#NJ)#XCB(NC 

3:C 


4-: 

DO 100 IC«1’#NC 

5: 

DO 200 I»1#NB 

6: 200 

XCB( IC# I)«0,0 

7: 

DO 300 I»1#NI 

8 : 

DO 300 J«1#NJ' 

9: 

ilB-IBd# J) 

lO: 

IF(IIB.EQ.O) GOTO 300 

li: 

XCBdC# IIB)«XCBdC» IIB)+R{‘I)*Z ( J# IC)+Bd# J) 

12: 300 

CONTINUE 

13: 100 

CONTINUE 

14: 

RETURN 

15 : 

END 


A6EN 


X: SUBROUTINE AGEN( A> NflA^NA^ NOA# NQl# IDQ> ZL#CSrSN» lU) 

2: DIMENSION A ( 1) / NMA'{ 2> U » TDO ( 1 ) # AS { 2# 2 ) # D { 2, 4 ) 

3:C 

4; DO 100 I«l>2 

3i 00 100 J*1^4 

6: D(I#J)-0.0 

7: DO 100 K»l,2 

bJ 100 AB(I^J»K)»0.0 
9: SDTO(10>20#30,40»50>60^70)* lU 

10:C 

11: 10 NA»2 

12: NQl-2 

13: NQ1B»1 

14: IDQ(i)»l 

15: IDQ(2>«2 

16: ABU>1>1)«1.0 

17: AB(1^2>1)— 1.0/ZL 

18: AB(2#2> 1)«1.0/ZL 

19: D(1#1)«CS 

20: 0(1>2)»SN 

21 : GOTO 200 

22: 20 NA»4 

23*. NQ1«3 

24: NQ1B«2 

25: IDQ(l)-3 

26: IDQ(2)»4 

27: I0Q(3)«5 


K2D-72 



28 

29 

30 

31 

32 

33 
3A 

35 

36 

37 

38 

39 
^0 

41 

42 

43 

44 

45 

46 

47 

48 

49 

50 

51 

52 

53 

54 

55 

56 

57 
56 

59 

60 
61 
62 

63 

64 

65 

66 
67 


ZL2-ZL+ZI 
ZL3-ZL2+ZL 
AB(1,1>U»1.0 
AB(1,3,D— 3.0/ZL2 
ABdf 4»1}»2.0/^L3 
A8(l#2,2)-1.0 
AB(l,3,2)»-2.0/Z{. 
AB(i^4#2)«1.0/ZL2 
AB<2»3 p 1)— AB<1# 3> 1) 
A8{2#4,1)—AB(1#4,1) 
AB(2»3,2)«-1.0/ZL 
AB(2#4>2)-1.0/ZL2 
D(l#l)»1.0 
D(2,2)«SN 
D(2^3)»^CS 
GOTO 200 
30 NA-2 
NQ1»2 
f401B>l 
ID0(1)«1 
iOQ(2)«2 
ABdf 1,1)«1.0 
AB(1,2>1)»-1,0/Zl 
AB(2/2d) »1,0/ZL 
D(1,1)-SN 
DU/2)— CS 
GOTO 200 
40 NA-2 
NQle2 
NQ1B»1 
I00(l)»4 
IDCt2)»5 
AB{ld#l)»1.0 
aB{ 1#2*1>— 1.0/ZL 
AB(2,2>1)»1.0/ZL 
DU,1)-CS 
D(1#2)»SN 
GOTO 200 
50 NA-1 
NOl-2 


68 5 
695 
70 : 
7i: 
72 5 
735 
74; 
755 


NQ1B«1 
IDCCl) -1 
XDQ(2) -2 

ABa»id)»1.0/ZL 
AB<2il#l»»-1.0/ZL 
D{l#i)»SN 
D(l»2) *-CS 
GOTO 200 


76: 60 NA»3 

77: NQl-3 

78: NQ1B«2 

795 IDQ(1)«3 

80: IDQ(2)»4 

bi: IDQ(3)»5 

82J ZL2«Z_L»^ 


K2D-73 



03: 


ZL3-ZL2+ZL 

84: 


AB(1>2^ 1)»“6,0/ZL2 

8bi 


AB{1,3^ D-6.0/ZL3 

86: 


AB(1>1, 2)»1.0 

87: 


A8(li2F2)»-4,0/ZL 

8o ! 


AB(1,3#2)»3,0/ZL2 

89: 


AB(2,2>1)«6.0/ZL2 

90: 


AB(2^ 3,1)— 6.0/ZL3 

9i: 


AB(2,2,-2)— 2.0/ZL 

92: 


AB(2,3,2)»3.0/ZL2 

93: 


D(l,l)-1.0 

94: 


D(2,2)«SN 

96: 


D(2,3)— CS 

96: 


GOTO 200 

97: 

70 

NA«2 

98: 


NQ1«1 

99: 


NQlB-1 

100 : 


IDQ( 1) «3 

iOi: 


AB(1,1,1)»1*0 

102: 


AB(1,2,1) — 1.0/ZL 

103: 


AB(2,2,1)»1.0/ZL 

104: 


D(l,l)»1.0 

105:C 



106: 

200 

CONTINUE 


107: NM»-1 

108: DO 300 I«1,NA 

109: Nh-Nf*+1 

110: nMA(1»I)«NM 

111: 300 N«A(2#I)«0 

112: NQA-2+NQ1 

113: CALL AID { A^ AB> 0, NA^ NCA# NQl B# NQl ) 

il^: RETURN 

li;j: END 


- ATD 


SUBROUTINE ATD ( A, A B, D,N A» NQA, NQIB # NQl ) 

DIMENSION A<NA«NOA)« AB ( 2> NQIS ) j> 0 ( 2> ^ ) 

JST— NQl 

DO 100 IEND«1,2 

JST«JST+NQ1 

DO 200 IA«1>NA 

DO 200 J»1>NQ1 

JJ-JST+J 

A(1A> JJ)>0.0 

DO 200 K»1,N01B 

A(IA# JJ)«A(IA> JJ)+AB{IEND^ IA>K)+0(K^ J) 
200 CONTINUE 
100 CONTINUE 
RETURN 
END 


2 

3 

‘i 

5 

6 
7 
6 
9 

10 

11 

12 

13 

14 
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NCALNA 


li SUBROUTINE NCALNA ( XL> NMC> NMA, NC, NA, NAL FA> AL F A, S IGN ) 

Zi DIMENSION Z(NC/NA)# NHC<2/NC)^Nf1A{2#NA)/ALFA(i) 

3: C0HM0N/CK2D5/IQX(10)>QX(10) 

4: DO 100 I-l/NC 

5: DO 100 J«1#NA 

6: N»NMC{l>n+NMA(l# J) 

7: IFCNALFA.GT.O) GOTO 200 

bi Nl-N+1 

9: IF(IQX(N1).6T*0) GOTO 400 


10; IQX(N1)»1 

11! AN1»N1 

12: QX(N1)»(XL**N1)/AN1 

13; 400 Z( J) »QX{N1) 

14; GOTO 100 

li»; 200 ZCIiJ)»0.0 

16; N5N-I 

17! DO 300 lAL-l/NALFA 

18; N«N+1 

19; N1«N+1 

20; IF( IQX(Nl) .GT.O) GOTO 300 

21: IQX(N1)»1 

22: AN1»N1 

23; QX(Nl) «{XL«*N1) /ANl 

24; 300 Z(I#J)«Z{I#J)+QX(N1)*ALFA(IAL)*SIGN 

25; 100 CONTINUE 

26; RETURN 

27: £N0 


CXA 


l: SUBROUTINE CXA ( C> X# A,T^ TMP# NA> NuA# NOT# NC# NBC# NBT# I NC B# MAPQ# I DQ 

2; *, N0DES#NIJ#NQ1> 

3; DIMENSION X(NC#NA)# A(NA#NQA)#T(NBT#NQT)#TMP{NBC#NA) 

4: C(NC#NBC)#IDQ(NQ1)#NIJ(NODES)#MAPO(4#1) 

5:C 

6!C*+» THIS ROUTINE COMPUTES T«C ( TRANSPOSE ) *X+A 

7:C+*+ RESULT IS ACCUMULATED. T MUST BE INITIALIZED OUTSIDE 

8:C 

9! DO 100 I»1#NBC 

10: 00 100 J«1>NA 
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11 ; 
12 ! 
13: 
14;C 
lij; 
16; 
17; 
18; 
19; 
ZOi 
21 ; 
22 ; 
23: 
24; 
25: 
26; 
27; 
28 ; 
29: 


TMPd, J)-0.0 
DO 100 K»1,NC 

100 TMPd, JKClKjfI)*X{K> J) 

00 200 I«1jNBC 
J*0 

IT-I+INCB 

DO 200 MD«l, NODES 

NODE»NIJ(ND) 

DO 200 IQ»l,NQl 
I1Q«ID0(1Q) 

J“J+1 

JT»KAPQ{NODb, IIQ) 

IFUT.LE.O) GO TO 200 
DO 300 K»1>NA 

300 TdT, JT)«T(IT^ JT) + TMP(IjK)*A(K> J) 
200 CONTINUE 
RETURN 
END 


. SSTM 


l: 
2 ; 
3; 
4 : 
5; 
bi 
7: 
: 

9: 

10: 

11 : 

12 : 

13: 

14: 

15; 

16: 


SUBROUTINE SSTM ( T> S# STM> NHD> NB^ NQ ) 

DIMENSION H(NHD,NHD),T(NHD>NO),Sd)/STh{NB,NQ) 
00 100 I»1#N8 
DO 100 J*1#NQ 
STMd^ J)«0. 0 
00 100 K«1»NB 

100 STMd, J)»STM(Id )+Hd;K)*T(K^ J) 

K-0 

00 200 I«1#NQ 
00 200 J»l»I 
K-K+1 
S(K)«0.0 
DO 200 L«1»NB 

200 S{K)»S{K) + T(Ld)+STMa/J) 

RETURN 

END 
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Expepzmentrl Element Crprpility 


Tp»LE DF CnNTENTE 


Section Contents' 

1.0 TERMINOLDGY 

2.0 USER tilRITTEN SUPPDUTINES 

3.0 RRD6PRM EXECUT];pN 


EXPE-I 



EXPE 1 . 0 


TEPMINDLD6Y 


FdP ERCH element » AN "ELEMENT PEFEPENCE FPRME" IS EEFINEH 
USING THE SRME CDNI/ENTIDN RS FDP STRNJJRPrr SPflP ELEMENT 

TYPES E3^* E41* S41» S61* rnb S81 » Thrt is» the x rxis is 

DIPECTED FPDM NODE 1 THPDUGH NODE £» RNI» NODE 3 LIES IN THE 
FI PST PURDPi^NT DF THE X“ Y PLRNE. 

The fdllpwing symbols will be useb in explrining how to employ 

THE EXPEPIMENTRL ELEMENT CRPREILITYl 

Hrme Dimension Definition 

Type R 4-chrprctep rlphrnumepic wdpe rssignei> 

BY THE USEP TO IDENTIFY R PRPTICULRP ELEMENT 
FDPMULRTIDNj RNRLOGDUS TO THE STRNDRPD 
ELEMENT TYPE CODES ^E£l y E43» S81 > ETC. > . 
fl MODEL MRY CONSIST OF RNY DF THE STRNDRPD 
ELEMENT TYPES? PLUS ONE OP MOPE TYPES OF 
EV-PEPIMENTRL ELEMENTS. ThEPE IS NO SPECIFIC 
LIMIT ON THE NUMBEP OF TYPES DF EXPEPIMENTRL 
ELEMENTS. 

MRJDP Elements rpe clrssified rs: 

MRJOP = 1 FDP LINE ELEMENTS* E.G. E£1 * E££. 

Mf=|JDP = £ FDP £D ELEMENTS* E.G. E31* E43 . 
MfiJDP E 3 FDP 3D ELEMENTS? E.G. 381? F41. 

MINOR Rn INTEGEP RESIGNED BY THE USEP TO IDENTIFY 

R SPECIFIC ELEMENT TYPE. R UNIOUE "mINDP" 

IS DEFINED FDP ERCH “TYPE. *' 

N The NUMBEP DF NODES PEP ELEMENT ^M IN I MUM = 3* 

MR VI mum = 3£> . EXPEPIMENTRL BERM ELEMENTS 
should be modelled rs 3-ndde elements* using 

THIPD POINT TO ESTRBLISH THE CROSS-SECTION 
DPIENTRTION. 

M N^N + l>x£. See K« CM* KG below. 

NDF MuMEEP DF DEGPEES DF FREEDOM PEP NODE? EITHER 

3 '' 3 DISPLRCEMENTS> * OP 6 < 3 DISPLRCEMENTS 

RND 3 PDTRTIDNS> . 

NNDF The TOTRL number of degrees of freedom pep 

ELEMENT? N TIMES NDF. SeE U RND F BELOW. 

X C3?N) PeLRTIUE TO the element reference FPRME* THE 

DIRECTION I POSITION CDDPDINRTE DF ELEMENT 
NODE J IS X <I ? . 
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TePMJNDLOiSY* CONTINUE!}! 

Nome DiMENffioff Iiefinitidn 


U <NN»F> 


F <HNr>F> 

P <’NP> 
n4ppdp 


S <NS) 


0 <NC>> 


Element nocf^l motion i/ectdp. Where i>za mni} 

f*lJ nPE riPECTIDN I I}I£PLfiCEMENT RND PDTftTIDN 
COMPONENTS <REURTIVE TO THE ELEMENT PEFEPENCE 
PPRME> DF ELEflENT NOPE J» THE OPPER OF TERMS 
IN U IS RS FOLLOWS? IF NPF=6S 

Dll p21 d 31 Rll PHI r31 p12 p22 d32 p12 

If NpF = 3> THE DPPEP IS THE SRME> EXCEPT 
THE PDTRTIDN TERMS RPE NOT PRESENT. 

Hone of the user— written subroutines rccess 
EITHER U OR F> PEFINED BELOW. U RND F RPE 
DEFINED HERE ONLY FDR USE IN DEFINING OTHER 
fiPRRYS WHICH MUST BE GENEPRTED IN USER- 
WRITTEN ROUTINES. 

Element nodrl force uectdr? corresponding 
TO U. 

Element property rprry. The content of P 
IS ESTRPLISHED BY THE USER. TyPICRL ITEMS 
RPE MRTERIRL RND SECTION PROPERTIES? OPTION 
CONTROLLERS? ETC. BEFORE EXECUTING ELD? 

THE USER MUST CPERTE? Ulfl flUS/TFiBLE? ft TftBUE 
NftMED XXXX BTRB 2 n4ppop? where n4ppdp 

IS RNY INTEGER GPEfiTER THRN 100. EftCH LINE 
IN THIS TftBLE IS ft P ftPPRY? ftPPLICRBLE TO 
ONE DP MORE INDIUIDUftL ELEMENTS? RS ESTAB- 
LISHED Ulfl THE McECT POINTER IN ELD. 

The ELEMENT "STRESS*' STATE VECTOR. The 
CONTENT OF S IS ESTABLISHED -BY THE USER. 
Typical terms in S -are stress field coeffi- 
cients? STRAINS? STRESS RESULTANTS? ETC. 

See the definition of O? P? and D below. 

The element "thermal" lord vector. The 
content of 0 IS established by the user. 
Typical components of 0 are coefficients 

DF TEMPERATURE DR DISLOCATION FUNCTIONS? 
TEMPERATURE GRADIENTS? ETC. ? DEFINED IN ANY 
MANNER THE USER FINDS CONVENIENT. SeE THE 
definition of the C and D MATRICES BELOW. 

The 0 VECTORS for individual elements are 
constructed VIA flUS/ELDfiTFI. The NODFIu 
temperature and HODRl PRESsure arrays ? if 

PRESENT* DO NOT APPLY TO EXPERIMENTAL 
ELEMENTS. 
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Tepmindlqgyj continued: 

MftME Dihension Definition 

Y ^N> Element nodal w^tisht dietpidutidn. The weight 

OF AN ELEMENT =*= Y <1> + V + “ - - Y <N> . 

I PpDCESSDP E USEE Y IN CONSTPUCTINe THE SYSTEM 

DIAGONAL MASS MATRIX? DEM. 

K fNDF?NDF»M> ThE ELEMENT STIFFNESS MATRIX? ARRAYED AS M 

SUPMATRZCES? EACH DIMENSIONED NDF PY NDF: 

K<-?-?l> KC-?-?£) KC-?-?4> - - - 

Kf~.-?3> K‘'-“*-?5> - - ^ 

K = K<-?-?6> - - - 

SYMMETRIC K<‘— ?-?m‘J 


C <NNDr?NP> Thermal force influence matrix. Total 

element nodal forces are: 

F = K U + CO 

P <ns?nndf'J Mechanical stress pecoi/epy matrix. 

D <’ns?np> Thermal stress recovery matrix. The Total 

stress? as computed py GSF? is: 

S - P U + DO 

CM <NDF?NDF»M> ELEMENT CONSISTENT MASS MATRIX? ARRAYED THE 

SAME AS K. 

KG Cndf?ndf?m‘> Element initial-stress stiffness matrix? 

ARRAYED THE SAME AS K. 
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EXPE £. 0 “ USER WRITTEN SUBROUTINES 


Tn IMPLEMENT EVREPlM^NTfiL ELEMENTEsi THE USER MUET CDZiE THE 
EOLLDMING ROUTINES HNr> INCPPPDPflTE THEM INTO THE INDICRTED 
PPDCESSDPS» PEPURCING EMPTY ROUTINES HRUING THE SRME NRMES 
IN THE STRNDRPP UEPSION OF SRftP. IT IS PERMISSIBLE TO OMIT 

DMEXPE> CMEXPEj rnd KGEXPE. 


Routine 

PpDCESSDP NrME fipGUMENT LIST 

inputs! outputs! 

E DMEXPE <■ MRJDP>MINDP«Nf NPF*P>X> Y> 

EKS KEXPE^ MR.TDP»MINDP*N*NPF<>P«X» K*PfC»B> 

M CMEXPE< mrjdpj minor* N«NDF'>P»X» CM^ 

KG KGEXPE < MRJORj MINOR? N* NPFyP*XjS* KG) 


The primary limitation on the eypepimentrl element capability 
IS the total length of an element information packet in the 
E“^tate> which is initiated in processor E. The component 

PARTS ARE! 

Part Length 

1 IP + N 

a 0 

3 9 + 1£n 

4 CP) NP 

5 ^'K) M X NDF X NDF 

6 CR) NS X NNDF 

7 (S) NS 

6 <C> NNDF X NP 

9 <D) NS X NP 

THE RELATED LIMITATIONS ARE! 

— The sum of the lengths of parts 1 through 4 may not 

EXCEED 1000. 

— fiLL PROCESSORS ACCESSING THE £~STATE <E.G.. K * M* KG» EONF * 
6SF) PERU I PE SUFFICIENT CENTRAL MEMORY TO LOAD AT LEAST 
ONE COMPLETE ELEMENT INFORMATION PACKET. • IN ADDITION TO 
THE OTHER C.M. REPUIPEMENTS? AS DEFINED IN THE SPfiP 

Peference Manual. 

— In E* THE ENTIRE XX^X BTRB £ n4ppdp array must be in 

CENTRAL MEMORY. 

RlTHDUGH the experimental element capability MRS BEEN STRUCTURED 
SUCH THAT THERE ARE FEW SIZE PESTPICTIDNSf USERS SHOULD BE 
AWARE THAT EXTREMELY LARGE ELEMENT INFORMATION PACKETS WILL 
GENERALLY RESULT IN PELATIUELY HIGH I/D COSTS. 
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EXPE 3.0 - PPOSPflM EXECUTION 


EyeCUTIDNS PPDCEEK the SRME PiS in R STRNJJftPP RNRLYSIS? EXCEPT 
THRT BEFORE EXECUTING ELU» THE UEEP MUST EXECUTE PUS/TflBLE TO 
CDNSTPUCT? FDP ERCH TYPE DF EXPEPIMENTRL Et-EMENT* R TRBLE 

NRHEP XXXX BTRP £ n4ppdp. The identifying integers n4ppop» 

MRY BE RNY NUMBER GPERTEP TMRN 100. ErCM CINE IN THIS TRBLE 
IS R P RPPRYj R£ PPEUIDUSLY DEFINED? CDNTRINING THE MRTEPIRL 
RND SECTION CDNSTRNTS RPPUICRBLE TO ONE OP MOPE ELEMENTS. 

It is rusd mrndrtopy thrt r MflTEPIflL CONSTANT drtr set be 

GENEPRTED IN TAB? Et-'EN THOUGH THE EXPEPIMENTRL ELEMENTS DO 
NOT REFERENCE NATC. ThIS IS PEPUIPED BECRUSE PROCESSOR E 

rutdmrticrlly lords NATC» rnd petrins it in centprl memory 

THROUGHOUT EXECUTION. 

ExpEPIMENTRL elements are defined in ELB RS ILLUSTRATED BELOH. 
Within the same ELB execution* other elehent types? either 

STANDARD OR EXPEPIMENTRL? MRY ALSO BE DEFINED. 

S»X0T ELD 

s. 

?• Define all elements of r giuen Types 

S 

EXPE Type? major? minor? n? ndf? ns? no* n4ppdpS 
r^SECT=K S Line k of XXXX BTAB 2 n4prdp is "P" fop 
1 elements defined subsepuently. 

j1 j£ - - jNl. First element connects joints j1 j£ - -jn. 
j1 j 2 - — JN’i Second element 


r^SECT-J 

S 

Jl j2 - 
j1 j2 - 


NSECT IS THE ONLY TRBLE POINTER COMMRND RPRLICRBLE TO EXPEPIMENTRL 
ELEMENTS. ThE I*10B* INC? RND GPOUP COMMANDS FUNCTION THE SRME RS 
FOR THE STANDARD ELEMENTS. ThE STANDARD MESH GENERATION FACIL- 
ITIES MRY BE USED FDP EXPEPIMENTRL ELEMENTS HRt^ING THE SRME 
NODAL ARRANGEMENT AS STANDARD ELEMENTS. 

If THERMAL LORDING IS PRESENT? AN ELD ATA— FORMAT DATA SET NAMED 

TEMP Type iset icase must be created via AUS/ELDATA before 
EONF IS executed: 

3X0T AUS 

ELDATA! TEMP Tyre iset icase 
S=1 

E=l: ol p2 — — ONr$ P<no items> for element 1? group 1. 
e=£: ol r2 ~ pnoS 0<no items) fop element £» group 1. 


$ Line j of XXXX BTAB 2 n4ppop is "P" for 

elements defined SUBSErUENTL"Y. 

- JN?. 

- JN?. 
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The functions pepfdpme» py SPflP ppdccesdps in euppdpt of 

EYFEPIHENTRL EUEMENTS RPE SUMMRPIZEH EELDW. MdnE DF THE DTHEP 
SPfiP FPDCEESOPS <E.6. THE FL.DT PPD<5PRH5P PECDGNI2E EXPEPIHENTfiL 
ELEMENTS. 

PpDCESSDP Function 

ELD PPDCESSES ELEMENT DEFINITIONS* RND CPEPTES DEF* GD» 

GTIT« PND DIP DRTR SETS. 

TDPD Accounts fop presence of expepimentpl elements in 

CDNSTPUCTING THE TDPOLDGICfiL PPPPYS* KMRP» AMfiP* 
ENftBLING SYSTEM MfiTPIX RSSEMBLY RNB FRCTOPINE TO 
BE PERFORMED. 

E Cpertes the F— strte dhtb sets fdp expepimentbl 

ELEMENTS* IN SKELETRL FORM. AdDS CONTRIBUTION 
OF EXPEPIMENTRL ELEMENT WEIGHTS TO DEM. 

EKS Inserts the K * P« C* rni> D rprrys into the E-strte. 

K Includes expepimentrl element K's in system K. 

M Includes expepimentrl element CM's in system CEM. 

KG Includes expepimentrl element KG''s in system KG. 

flUt Permits lording of TEMP Type drtr sets i^ir ELDRTR. 

EC'MF Adds contribution of expepimentrl element fixed- 

joint FORCES TO ERUIURLENT NODRL FORCE RPRRYS. 

GoF Computes S for erch element* including thepmrl 

EFFECTS* IF PRESENT. S IS EMEEDDED IN THE E— STRTE* 
FDR USE BY KG. IF REPUESTED l•>Ifi THE EMBED RESET 
CONTROL- S IS NOT PRINTED BY GS'F UNLESS OMLIME=c‘. 

CTPS Type drtr set rpe not produced by GSF for 

EXPEPIMENTRL ELEMENTS. P?F DDES NOT RECOGNIZE 
EVPEPIMENTRL ELEMENTS. 
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The functidne pepfdpmeij py SPflP ppdceesdpe in support of 

EYPEPlHENTftL EL.EMENTE PRC SUMMARIZED EFL.DH. HdNE DF THE OTHER 
SPflP PROCESSORS <E.G. THE PLOT PRDCPAMS^ RECOGNIZE EXPER IMENTAL 
ELEMENTS. 

Processor Function 

ELI' Processes element definitions* pnp creates DEFj GD> 

GTIT* RND DIP DATA SETS. 

TQPD fiCCDUNTS FDR PRESENCE DF EXPERIMENTAL ELEMENTS IN 

CONSTRUCTING THE TOPOLOGICAL ARRAYS* KMflP* flPIfiP" 
ENABLING SYSTEM MATRIX ASSEMBLY ANP FACTORING TO 
BE PERFORMED. 

E Creates the F-state data sets for experimental 

ELEMENTS* IN SKELETAL FORM. RDDS CONTRIBUTION 
DF EXPERIMENTAL ELEMENT WEIGHTS TO DEM. 

EK? If.|SERTS THE K» P* C> AND D ARRAYS INTO THE E~STATE. 

K Includes experimental element K's in system K. 

M Includes experimental element CM''s in system CEM. 

KG Includes experimental element KG^s in system KG. 

flUS Permits loading of TEMP Type data sets ••'ia ELDflTfl. 

EONF Adds contribution of experimental element fixed- 

joint FORCES TO ERUtUALENT NODAL FORCE ARRAYS. 

GSF Computes S fop each element ' including thermal 

EFFECTS* IF PRESENT. S' IS EMBEDDED IN THE E~STATE' 
FDR USE BY KG* IF REQUESTED I'lA THE EMBED RESET 
CONTROL. S IS NOT PRINTED BY GSF UNLESS DNLINE=c:. 

STPS Type data set are not produced by GSF fop 

EXPERIMENTAL ELEMENTS . PSF DOES NOT RECOGNIZE 
experimental ELEMENTS. 
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